Cavity-induced emergent topological spin textures in a Bose–Einstein condensate

The experimental progresses in reaching the quantum degeneracy limit in atomic gases paved the way for the realization of quantum many-body phenomena in these highly tunable systems [1, 2]. Some of the most remarkable examples include the realization of the superfluid to Mott-insulator quantum phase transition [3], quantum magnetism and magnetic orderings [4–8], synthetic magnetic fields (i.e. Abelian gauge potentials) [9–11], and spin–orbit coupling (i.e. non-Abelian gauge potentails) [12–14] in ultracold quantum gases. While the first generation of experiments was limited to static lattices and local contact interactions, the study of highly nonlinear optical systems, where the back-action of the matter on the radiation fields is not negligible, has opened up new frontiers towards dynamical optical potentials, long-range atom-atom interactions, and exotic collective phenomena [15–19]. The most prominent examples include the coupling of ultracold atoms to dynamic quantized radiation fields of high-quality cavities [20–22], leading to the realization of the Dicke superradiance phase transition [23, 24], atomic recoil lasing [25–27], and the quantum phase transition between superfluid, supersolid, Mott-insulator, and density-wave phases [28, 29].

Almost all experimental works and most theoretical studies of coupled atom-cavity systems in the past were limited to systems where either the atomic internal states [30–33] or the atomic external degrees of freedom [34–45] are taken into account. Only recently theoretical investigations towards including both atomic internal and external degrees of freedom in cavity QED have been conducted. These systems exhibit many more interesting phenomena, including the emergence of synthetic strong magnetic fields and spin–orbit coupling [46–49], disorder-driven density and spin self-ordering [50], topological states [51, 52], and a variety of magnetic orders [53, 54]. Very recently, first experimental implementations showing spin self-ordering in a spin-1 BEC in a standing wave cavity were realized [55, 56].

In this work we study a novel type of hybrid atom-cavity system, where Λ-type spin-1 ultracold bosons are confined into quasi one dimension along one leg of a ring cavity with two 'pairs' of nearly resonant modes [57, 58]. Each pair consists of two counterpropagating modes with the same polarization which is orthogonal to the polarization of the other pair. In contrast to [59], each pair of modes only couples to one of the transitions in the Λ scheme, where the atoms are assumed to posses two ground states (e.g. two different Zeeman sublevels) and an electronic excited state; see figure 1. The adiabatic elimination of the upper atomic electronic state as well as of other Zeeman sublevels results in an effective two-component pseudospin model. We consider a case where two counterpropagating modes of orthogonal polarizations are externally pumped by lasers through the cavity mirrors as depicted in figure 1. These two modes do not interfere with one another and the system is, therefore, initially homogeneous. This is reminiscent of the scheme for generating equal Rashba–Dresselhaus synthetic spin–orbit coupled for neutral atoms in free space [12].

Zoom In Zoom Out Reset image size

That said, in the present case the cavity modes are dynamic and affected by the atomic dynamics as well as by photon losses through the cavity mirrors. Crucially, the unpumped modes are 'dynamically' populated by coherently scattered photons and, therefore, couple to the pumped modes and the atomic internal and external degrees of freedom. As an important consequence, the cavity-induced spin–orbit coupling for the atoms only emerges above a critical pump power, which in turn gives rise to novel nonequilibrium quantum phases and quantum phase transitions of various natures in our system. It is this dynamical population of the unpumped cavity modes and its nontrivial interplay with the other degrees of freedom which marks a sharp contrast to the free-space spin–orbit coupling [12–14, 60, 61] as well as all other previous cavity-based spin–orbit coupling schemes [46–48, 51].

As the frequency and the strength of the pump lasers are varied, in the mean-field regime the system displays three fundamentally different phases with distinct characteristics in density and spatial spin texture as shown in figure 2. The first phase is the density-wave–spin-wave (DW-SW) state, where the density has a crystalline structure and the pseudospin exhibits a spatial spin-wave texture along the cavity axis; see figure 4(a). The second phase is the plane-wave–spin-spiral (PW-SS) state, where the density is homogeneous while the pseudospin exhibits a spin-spiral texture as illustrated in figure 4(b). The third phase is the density-wave–spin-spiral (DW-SS) state, where a crystalline-ordered density is accompanied with a spin-spiral pseudospin texture; see figure 4(c). The latter two states with the spin-spiral texture have a Skyrmionic nature with a nontrivial topology [62, 63], and are intimately related to the emergence of cavity-induced spin–orbit coupling. The topologically trivial DW-SW state, on the other hand, breaks the continuous screw-like symmetry of the system, resulting in the appearance of a gapless Goldstone mode in addition to the phonon sound mode, a characteristic of a supersolid with two broken continuous symmetries [14, 29, 45]. The topological phase transitions between the DW-SW and the spin-spiral states, which is a direct consequence of the emergence of cavity-induced spin–orbit coupling, exhibit first-order characteristics, while the quantum phase transition between the two spin-spiral states are second-order. Remarkably, all the quantum phase transitions can be monitored in situ through the cavity output, as can be seen from the inset of figure 2.

Zoom In Zoom Out Reset image size

The paper is organized as follows. We introduce the model in section 2 and then derive the effective Hamiltonian and the Heisenberg equations of motion. We then find the steady-state solutions of the equations of motion in the mean-field limit in section 3. In this section we discuss the atomic phase diagram (section 3.1), the cavity-field phase diagram (section 3.2), and the effect of cavity-induced emergent spin–orbit coupling (section 3.3). Section 4 is devoted to the elementary excitations and the (broken) symmetries of the system. We present the concluding remarks in section 5. Appendices A and B show the details of the adiabatic elimination of the atomic excited state and linearization of the Heisenberg equations of motion, respectively.

Consider a Λ-type spin-1 BEC tightly confined into quasi-one dimension along one leg of a ring cavity in z direction as depicted in figure 1. The internal atomic states of interest consist of two pseudospin ground states, designated by $\left.\left|\downarrow \right.,\rangle \right\rangle$ and $\left.\left|\uparrow \right.,\rangle \right\rangle$, and an electronic excited state $\left|e\right\rangle$ with energies $\left\{{ \hbar }{\omega }_{\downarrow },{ \hbar }{\omega }_{\uparrow },{ \hbar }{\omega }_{e}\right\}$. The transition $\left.\left|\downarrow \right.,\rangle \right\rangle \leftrightarrow \left|e\right\rangle$ ($\left.\left|\uparrow \right.,\rangle \right\rangle \leftrightarrow \left|e\right\rangle$) couples to a pair of degenerate, counterpropagating electromagnetic modes ${\hat{a}}_{\pm }{{\rm{e}}}^{\pm {{\rm{i}}{k}}_{a}z}$ (${\hat{b}}_{\pm }{{\rm{e}}}^{\pm {{\rm{i}}{k}}_{b}z})$ of the ring cavity as shown in figure 1. The operator ${\hat{a}}_{+/-}$ (${\hat{b}}_{+/-}$) annihilates a forward/backward moving photon in the first (second) pair of ring-cavity modes with a wave vector ${k}_{a}={\omega }_{a}/c=2\pi /{\lambda }_{a}$ (${k}_{b}={\omega }_{b}/c=2\pi /{\lambda }_{b}$). For our desired system, the condition $| {\omega }_{a}-{\omega }_{b}| /{\omega }_{a(b)}\ll 1$ holds in general. Therefore, the assumptions $k:= {k}_{a}\approx {k}_{b}$ and $\lambda := {\lambda }_{a}\,\approx {\lambda }_{b}$ are legitimate and will be used throughout this work. Each pair of modes (i.e. ${\hat{a}}_{\pm }$ or ${\hat{b}}_{\pm }$) has the same polarization, which is orthogonal to the polarization of the other pair; e.g. $\left\{{\epsilon }_{+},{\epsilon }_{-}\right\}$. The mode ${\hat{a}}_{+}$ (${\hat{b}}_{-}$) is driven by an external pump laser through the cavity mirror with a frequency ${\tilde{\omega }}_{a}$ (${\tilde{\omega }}_{b}$) and an amplitude ${\eta }_{+}$ (${\eta }_{-}$).

The single-particle Hamiltonian density of the system in the dipole and rotating-wave approximation reads

$$\begin{eqnarray}&&{ \mathcal H }={{ \mathcal H }}_{0,\mathrm{at}}+{H}_{0,\mathrm{cav}}+{{ \mathcal H }}_{\mathrm{ac}}+{H}_{\mathrm{pump}},\end{eqnarray} \tag{ 1 }$$

with

$$\begin{eqnarray}&&{{ \mathcal H }}_{0,\mathrm{at}}=\displaystyle \frac{{\hat{p}}^{2}}{2m}+{ \hbar }\displaystyle \sum _{\tau =\downarrow ,\uparrow ,e}{\omega }_{\tau }{\hat{\sigma }}_{\tau \tau },\end{eqnarray} \tag{ 2a }$$

$$\begin{eqnarray}&&{H}_{0,\mathrm{cav}}={ \hbar }\displaystyle \sum _{j=+,-}\left({\omega }_{a}{\hat{a}}_{j}^{\dagger }{\hat{a}}_{j}+{\omega }_{b}{\hat{b}}_{j}^{\dagger }{\hat{b}}_{j}\right),\end{eqnarray} \tag{ 2b }$$

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{ac}}={ \hbar }\left\{\left[{{ \mathcal G }}_{\downarrow }\left({{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{a}}_{+}+{{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{a}}_{-}\right){\hat{\sigma }}_{{\rm{e}}\downarrow }+{{ \mathcal G }}_{\uparrow }\left({{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{b}}_{+}+{{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{b}}_{-}\right){\hat{\sigma }}_{{e}\uparrow }\right]+{\rm{h}}.{\rm{c}}.\right\},\end{eqnarray} \tag{ 2c }$$

$$\begin{eqnarray}&&{H}_{\mathrm{pump}}={\rm{i}}{ \hbar }\left({\eta }_{+}{\hat{a}}_{+}^{\dagger }{{\rm{e}}}^{-{\rm{i}}{\tilde{\omega }}_{a}t}+{\eta }_{-}{\hat{b}}_{-}^{\dagger }{{\rm{e}}}^{-{\rm{i}}{\tilde{\omega }}_{b}t}-{\rm{h}}.{\rm{c}}.\right).\end{eqnarray} \tag{ 2d }$$

Here m is the atomic mass, $\hat{p}=-{\rm{i}}\hslash {{\rm{\partial }}}_{z}$ is the atomic momentum operator and ${\hat{\sigma }}_{\tau \tau ^{\prime} }=\left|\tau \right\rangle \left\langle \tau ^{\prime} \right|$ are the atomic transition operators. The atom-photon coupling for the transition $\left.\left|\downarrow \right.,\rangle \right\rangle \leftrightarrow \left|e\right\rangle$ ($\left.\left|\uparrow \right.,\rangle \right\rangle \leftrightarrow \left|e\right\rangle$) is denoted as ${{ \mathcal G }}_{\downarrow }$ (${{ \mathcal G }}_{\uparrow }$) and h.c. stands for the Hermitian conjugate. ${{ \mathcal H }}_{0,\mathrm{at}}$ and ${H}_{0,\mathrm{cav}}$ are the bare atomic Hamiltonian density and the cavity-field Hamiltonian, respectively. ${{ \mathcal H }}_{\mathrm{ac}}$ represents the coupling between the atom and the cavity fields, and ${H}_{\mathrm{pump}}$ accounts for the pumping of the cavity fields.

The corresponding many-body Hamiltonian is obtained as

$$\begin{eqnarray}&&H=\int {\hat{{\rm{\Psi }}}}^{\dagger }(z)({{ \mathcal H }}_{0,\mathrm{at}}+{{ \mathcal H }}_{\mathrm{ac}})\hat{{\rm{\Psi }}}(z){\rm{d}}{z}+{H}_{0,\mathrm{cav}}+{H}_{\mathrm{pump}}+{H}_{\mathrm{int}},\end{eqnarray} \tag{ 3 }$$

where $\hat{{\rm{\Psi }}}={\left({\hat{\psi }}_{e},{\hat{\psi }}_{\uparrow },{\hat{\psi }}_{\downarrow }\right)}^{T}$ is the three-component atomic field operator which satisfies the bosonic commutation relation $[{\hat{\psi }}_{\tau }(z),{\hat{\psi }}_{\tau ^{\prime} }^{\dagger }(z^{\prime} )]=\delta (z-z^{\prime} ){\delta }_{\tau ,\tau ^{\prime} }$. The Hamiltonian ${H}_{\mathrm{int}}$ accounts for two-body contact interactions between the atoms, and ensures the thermalization and relaxation of the BEC. However, we assume that the two-body contact interactions are negligibly small compared to cavity-mediated long-range interactions, which is a good approximation for typical cavity-QED experiments [29, 64]. Therefore, we will not explicitly include two-body contact interactions in our model.

The dynamics of the system is governed by the Heisenberg equations of motion of the atomic field operators ${\rm{i}}{ \hbar }{\partial }_{t}{\hat{\psi }}_{\tau }=[{\hat{\psi }}_{\tau },H]$, and the photonic field operators ${\rm{i}}{ \hbar }{\partial }_{t}{\hat{a}}_{\pm }=[{\hat{a}}_{\pm },H]-{\rm{i}}{ \hbar }\kappa {\hat{a}}_{\pm }$ and ${\rm{i}}{ \hbar }{\partial }_{t}{\hat{b}}_{\pm }=[{\hat{b}}_{\pm },H]-{\rm{i}}{ \hbar }\kappa {\hat{b}}_{\pm }$. The decay (i.e. leakage) of cavity photons is included phenomenologically by adding the terms proportional to κ in the latter equations for the photonic-field operators. We assume equal decay rates for the two modes with orthogonal polarizations which, in general, is not entirely true in experimental setups because the cavity mirrors are polarization selective. However, different decay rates only lead to quantitative changes in the results presented in the following, whereas the qualitative physics does not change. If the relative atomic detunings with respect to the pump lasers ${{\rm{\Delta }}}_{\downarrow (\uparrow )}:= {\tilde{\omega }}_{a(b)}-[{\omega }_{e}-{\omega }_{\downarrow (\uparrow )}]$ are large compared to the two-photon detuning $\delta := {{\rm{\Delta }}}_{\uparrow }-{{\rm{\Delta }}}_{\downarrow }$ and the atom-photon couplings $\left\{{{ \mathcal G }}_{\downarrow },{{ \mathcal G }}_{\uparrow }\right\}$, the atomic excited state reaches a steady-state on a short time scale and its dynamics can be adiabatically eliminated. This results in a set of six coupled effective Heisenberg equations for the atomic pseudospin and photonic field operators

$$\begin{eqnarray}&&{\rm{i}}{ \hbar }\displaystyle \frac{\partial }{\partial t}\left(\begin{array}{c}{\hat{\psi }}_{\downarrow }\\ {\hat{\psi }}_{\uparrow }\end{array}\right)={{ \mathcal H }}_{\mathrm{at}}\left(\begin{array}{c}{\hat{\psi }}_{\downarrow }\\ {\hat{\psi }}_{\uparrow }\end{array}\right),\quad \,{\rm{i}}{ \hbar }\displaystyle \frac{\partial }{\partial t}\left(\begin{array}{c}{\hat{a}}_{+}\\ {\hat{a}}_{-}\\ {\hat{b}}_{+}\\ {\hat{b}}_{-}\end{array}\right)={{ \mathcal H }}_{\mathrm{cav}}\left(\begin{array}{c}{\hat{a}}_{+}\\ {\hat{a}}_{-}\\ {\hat{b}}_{+}\\ {\hat{b}}_{-}\end{array}\right)+{\rm{i}}{ \hbar }\left(\begin{array}{c}{\eta }_{+}\\ 0\\ 0\\ {\eta }_{-}\end{array}\right).\end{eqnarray} \tag{ 4 }$$

The details are presented in appendix A.

Now ${{ \mathcal H }}_{\mathrm{at}}$ and ${{ \mathcal H }}_{\mathrm{cav}}$ are the 'effective' atomic and cavity-field Hamiltonian densities, respectively, which contain the couplings between all atomic and photonic degrees of freedom. The effective atomic Hamiltonian density has the matrix form

$$\begin{eqnarray}{{ \mathcal H }}_{\mathrm{at}}=\left(\begin{array}{cc}\tfrac{{\hat{p}}^{2}}{2m}+{ \hbar }{\hat{U}}_{\downarrow }(z)-\tfrac{{ \hbar }\delta }{2} & { \hbar }{\hat{{\rm{\Omega }}}}_{{\rm{R}}}(z)\\ { \hbar }{\hat{{\rm{\Omega }}}}_{{\rm{R}}}^{\dagger }(z) & \tfrac{{\hat{p}}^{2}}{2m}+{ \hbar }{\hat{U}}_{\uparrow }(z)+\tfrac{{ \hbar }\delta }{2}\end{array}\right),\end{eqnarray} \tag{ 5 }$$

with the potential operators

$$\begin{eqnarray}&&{\hat{U}}_{\downarrow }(z)={U}_{0\downarrow }\left({\hat{a}}_{+}^{\dagger }{\hat{a}}_{+}+{\hat{a}}_{-}^{\dagger }{\hat{a}}_{-}+{{\rm{e}}}^{-2{\rm{i}}{k}{z}}{\hat{a}}_{+}^{\dagger }{\hat{a}}_{-}+{{\rm{e}}}^{2{\rm{i}}{k}{z}}{\hat{a}}_{-}^{\dagger }{\hat{a}}_{+}\right),\end{eqnarray} \tag{ 6a }$$

$$\begin{eqnarray}&&{\hat{U}}_{\uparrow }(z)={U}_{0\uparrow }\left({\hat{b}}_{+}^{\dagger }{\hat{b}}_{+}+{\hat{b}}_{-}^{\dagger }{\hat{b}}_{-}+{{\rm{e}}}^{-2{\rm{i}}{k}{z}}{\hat{b}}_{+}^{\dagger }{\hat{b}}_{-}+{{\rm{e}}}^{2{\rm{i}}{k}{z}}{\hat{b}}_{-}^{\dagger }{\hat{b}}_{+}\right),\end{eqnarray} \tag{ 6b }$$

and the two-photon Raman coupling operator

$$\begin{eqnarray}&&{\hat{{\rm{\Omega }}}}_{{\rm{R}}}(z)={{\rm{\Omega }}}_{0{\rm{R}}}\left({\hat{a}}_{+}^{\dagger }{\hat{b}}_{+}+{\hat{a}}_{-}^{\dagger }{\hat{b}}_{-}+{{\rm{e}}}^{-2{\rm{i}}{k}{z}}{\hat{a}}_{+}^{\dagger }{\hat{b}}_{-}+{{\rm{e}}}^{2{\rm{i}}{k}{z}}{\hat{a}}_{-}^{\dagger }{\hat{b}}_{+}\right).\end{eqnarray} \tag{ 7 }$$

Here we have introduced the maximum depth of the optical potential per photon ${U}_{0\tau }:= 2| {{ \mathcal G }}_{\tau }{| }^{2}/({{\rm{\Delta }}}_{\downarrow }+{{\rm{\Delta }}}_{\uparrow })$ and the maximum two-photon Raman transition frequency ${{\rm{\Omega }}}_{0{\rm{R}}}:= 2{{ \mathcal G }}_{\downarrow }^{* }{{ \mathcal G }}_{\uparrow }/({{\rm{\Delta }}}_{\downarrow }+{{\rm{\Delta }}}_{\uparrow })$. The former potential depth results from two-photon scatterings between cavity modes with the same polarization without changing the atomic internal state, whereas the latter frequency ${{\rm{\Omega }}}_{0{\rm{R}}}$ is due to two-photon scatterings between cavity modes with orthogonal polarizations accompanied by an atomic pseudospin flip $\left.\left|\downarrow \right.,\rangle \right\rangle \leftrightarrow \left.\left|\uparrow \right.,\rangle \right\rangle$.

The coherent scattering of photons from a mode back into itself results in the position-independent energy shifts (i.e. the terms ${\hat{a}}_{\pm }^{\dagger }{\hat{a}}_{\pm }$ and ${\hat{b}}_{\pm }^{\dagger }{\hat{b}}_{\pm }$) in the potential operators ${\hat{U}}_{\tau }(z)$, while the scattering of photons from a mode (say, ${\hat{a}}_{+}$) into its degenerate counterpropagating mode (i.e. ${\hat{a}}_{-}$) results in the position-dependent terms proportional to ${{\rm{e}}}^{\pm 2{\rm{i}}{k}{z}}$ in the potential operators. The former photon scatterings do not transfer any momentum to the atom, while the latter scatterings lead to $\pm 2{ \hbar }k$ momentum kicks to the atom. These particular processes are schematically shown in figures 3(a) and (b), respectively.

Zoom In Zoom Out Reset image size

The photon scattering between modes with orthogonal polarizations and the same propagation direction (that is, between ${\hat{a}}_{\pm }$ and ${\hat{b}}_{\pm }$) gives rise to the position-independent Raman coupling terms (i.e. the terms ${\hat{a}}_{+}^{\dagger }{\hat{b}}_{+}$ and ${\hat{a}}_{-}^{\dagger }{\hat{b}}_{-}$) in ${\hat{{\rm{\Omega }}}}_{{\rm{R}}}(z)$, while scattering between modes with orthogonal polarizations and opposite propagation directions (that is, between ${\hat{a}}_{\pm }$ and ${\hat{b}}_{\mp }$) results in the position-dependent Raman coupling terms proportional to ${{\rm{e}}}^{\pm 2{\rm{i}}{k}{z}}$ in the two-photon Raman coupling operator. These two-photon Raman processes are illustrated in figures 3(c) and (d). While the atomic pseudospin is flipped in both Raman processes, there is no net momentum transfer to the atom in the former processes, whereas the latter processes impart a ±2${ \hbar }$k momentum to the atom. These latter Raman transition processes with ±2${ \hbar }$k momentum kicks can induce a synthetic spin–orbit coupling for the atom.

The effective cavity-field Hamiltonian 'density' in matrix form is given by

$$\begin{eqnarray}{{ \mathcal H }}_{\mathrm{cav}}={ \hbar }\left(\begin{array}{cccc}-{\tilde{{\rm{\Delta }}}}_{a} & {U}_{0\downarrow }{\hat{{ \mathcal N }}}_{\downarrow } & {{\rm{\Omega }}}_{0{\rm{R}}}{\hat{S}}_{-} & {{\rm{\Omega }}}_{0{\rm{R}}}{\hat{{ \mathcal S }}}_{-}^{(-)}\\ {U}_{0\downarrow }{\hat{{ \mathcal N }}}_{\downarrow }^{\dagger } & -{\tilde{{\rm{\Delta }}}}_{a} & {{\rm{\Omega }}}_{0{\rm{R}}}{\hat{{ \mathcal S }}}_{-}^{(+)} & {{\rm{\Omega }}}_{0{\rm{R}}}{\hat{S}}_{-}\\ {{\rm{\Omega }}}_{0{\rm{R}}}^{* }{\hat{S}}_{+} & {{\rm{\Omega }}}_{0{\rm{R}}}^{* }{\hat{{ \mathcal S }}}_{+}^{(-)} & -{\tilde{{\rm{\Delta }}}}_{b} & {U}_{0\uparrow }{\hat{{ \mathcal N }}}_{\uparrow }\\ {{\rm{\Omega }}}_{0{\rm{R}}}^{* }{\hat{{ \mathcal S }}}_{+}^{(+)} & {{\rm{\Omega }}}_{0{\rm{R}}}^{* }{\hat{S}}_{+} & {U}_{0\uparrow }{\hat{{ \mathcal N }}}_{\uparrow }^{\dagger } & -{\tilde{{\rm{\Delta }}}}_{b}\end{array}\right),\end{eqnarray} \tag{ 8 }$$

where we have introduced the effective cavity detunings ${\tilde{{\rm{\Delta }}}}_{a(b)}:= ({{\rm{\Delta }}}_{a(b)}+{\rm{i}}\kappa )-{U}_{0\downarrow (\uparrow )}{\hat{N}}_{\downarrow (\uparrow )}$ with the relative cavity detunings with respect to the pump frequencies ${{\rm{\Delta }}}_{a(b)}:= {\tilde{\omega }}_{a(b)}-{\omega }_{a(b)}$ and the particle number operators

$$\begin{eqnarray}&&{\hat{N}}_{\tau }=\int {\hat{\psi }}_{\tau }^{\dagger }(z){\hat{\psi }}_{\tau }(z){\rm{d}}{z}.\end{eqnarray} \tag{ 9 }$$

The off-diagonal coupling operators are given by

$$\begin{eqnarray}&&{\hat{{ \mathcal N }}}_{\tau }=\int {{\rm{e}}}^{-2{\rm{i}}{k}{z}}{\hat{\psi }}_{\tau }^{\dagger }(z){\hat{\psi }}_{\tau }(z){\rm{d}}{z},\end{eqnarray} \tag{ 10a }$$

$$\begin{eqnarray}&&{\hat{S}}_{-}={\hat{S}}_{+}^{\dagger }=\int {\hat{\psi }}_{\downarrow }^{\dagger }(z){\hat{\psi }}_{\uparrow }(z){\rm{d}}{z},\end{eqnarray} \tag{ 10b }$$

$$\begin{eqnarray}&&{\hat{{ \mathcal S }}}_{-}^{(\pm )}={\left({\hat{{ \mathcal S }}}_{+}^{(\mp )}\right)}^{\dagger }=\int {{\rm{e}}}^{\pm 2{\rm{i}}{k}{z}}{\hat{\psi }}_{\downarrow }^{\dagger }(z){\hat{\psi }}_{\uparrow }(z){\rm{d}}{z},\end{eqnarray} \tag{ 10c }$$

where ${\hat{S}}_{+}$ (${\hat{S}}_{-}$) is the collective atomic spin raising (lowering) operator. The operators ${\hat{{ \mathcal N }}}_{\tau }$ and ${\hat{{ \mathcal S }}}_{\pm }^{(\pm )}$ are the density- (for pseudospin τ) and spin-wave operators, respectively.

The matrix elements of the effective cavity-field Hamiltonian density ${{ \mathcal H }}_{\mathrm{cav}}$ give the strengths of the two-photon processes depicted in figure 3. The diagonal terms proportional to ${\hat{N}}_{\tau }$ correspond to the processes illustrated in figure 3(a) and result in the dispersive shifts ${U}_{0\tau }{\hat{N}}_{\tau }$ of the cavity frequencies. Terms proportional to ${\hat{{ \mathcal N }}}_{\tau }$ (and their Hermitian conjugates) correspond to the processes shown in figure 3(b) and provide the couplings between degenerate counterpropagating modes with the same polarization. Finally, the matrix elements proportional to ${\hat{S}}_{\pm }$ and ${\hat{{ \mathcal S }}}_{\pm }^{(\pm )}$ correspond to the pseudospin flipping processes depicted in figures 3(c) and (d), respectively, and yield the nontrivial couplings between modes with orthogonal polarizations. Note that the relative time phase between the two pump fields has to be fixed in order to drive the coherent processes described above.

The system possesses a screw-like continuous symmetry. It is manifested in the invariance of the total effective Hamiltonian corresponding to (4) under a simultaneous spatial translation $z\to z+{\rm{\Delta }}z$, phase rotations of the unpumped photonic field operators ${\hat{a}}_{-}\to {\hat{a}}_{-}{{\rm{e}}}^{2{\rm{i}}{k}{\rm{\Delta }}z}$ and ${\hat{b}}_{+}\to {\hat{b}}_{+}{{\rm{e}}}^{-2{\rm{i}}{k}{\rm{\Delta }}z}$, and phase rotations of the atomic field operators ${\hat{\psi }}_{\downarrow }\to {\hat{\psi }}_{\downarrow }{{\rm{e}}}^{-{\rm{i}}{k}{\rm{\Delta }}z}$ and ${\hat{\psi }}_{\uparrow }\to {\hat{\psi }}_{\uparrow }{{\rm{e}}}^{+{\rm{i}}{k}{\rm{\Delta }}z}$. Note that the phases of the pumped cavity modes ${\hat{a}}_{+}$ and ${\hat{b}}_{-}$ are fixed by the cavity pumps η_±. The phase rotation of ${\hat{a}}_{-}$ and ${\hat{b}}_{+}$ results in a shift of the potential minima in (6), defining the position of the atomic-density maxima. On the other hand, the phase shift of the atomic field operators leads to rotations of the spin operators: ${\hat{S}}_{+}\to {\hat{S}}_{+}{{\rm{e}}}^{-2{\rm{i}}{k}{\rm{\Delta }}z}$ and ${\hat{{ \mathcal S }}}_{+}^{(\pm )}\to {\hat{{ \mathcal S }}}_{+}^{(\pm )}{{\rm{e}}}^{-2{\rm{i}}{k}{\rm{\Delta }}z}$. Hence, a translation in space is tied to a corresponding rotation of the atomic spin, where the rotation angle directly depends on the length of the spatial translation Δz. This leads to a screw-like continuous symmetry.

The system is highly nonlinear. The effective atomic Hamiltonian density depends on the cavity fields through the potential ${\hat{U}}_{\tau }(z)$ and the Raman operators ${\hat{{\rm{\Omega }}}}_{{\rm{R}}}(z)$, while the effective cavity-field Hamiltonian density depends on the atomic fields via the atomic number ${\hat{N}}_{\tau }$, the density-wave ${\hat{{ \mathcal N }}}_{\tau }$, the collective atomic spin ${\hat{S}}_{\pm }$, and the spin-wave ${\hat{{ \mathcal S }}}_{\pm }^{(\pm )}$ operators. It is this nonlinear dynamics and the nontrivial interplay between various degrees of freedom which give rise to intriguing phenomena in our system, as it will be discussed in the subsequent sections. These nonlinear dynamics and the nontrivial interplay between the cavity modes survive even in the strong pumping limit ${\eta }_{\pm }\gg {\omega }_{\mathrm{rec}}$, with ${\omega }_{\mathrm{rec}}:= { \hbar }{k}^{2}/2m$ being the recoil frequency. Although in the strong pumping limit the pumped cavity fields $\left\{{\hat{a}}_{+},{\hat{b}}_{-}\right\}$ behave as classical fields, the unpumped modes $\left\{{\hat{a}}_{-},{\hat{b}}_{+}\right\}$ still retain their quantum nature and behave as dynamical fields.

In the following, we restrict our analysis to red-detuned pump lasers with respect to both bare atomic and cavity frequencies, i.e. $\left\{{{\rm{\Delta }}}_{\downarrow ,\uparrow },{{\rm{\Delta }}}_{a,b}\right\}\lt 0$. The atoms are therefore attracted to the intensity maxima of the light fields, while experiencing cavity cooling. In order to reduce the number of free parameters and capture the fundamental physics, we further focus on the special case of a completely symmetric configuration: ${\rm{\Delta }}:= {{\rm{\Delta }}}_{a}={{\rm{\Delta }}}_{b}$, $\eta := {\eta }_{+}={\eta }_{-}$, δ = 0, and ${{ \mathcal G }}_{\downarrow }={{ \mathcal G }}_{\uparrow }$ which results in ${U}_{0}:= {U}_{0\downarrow }={U}_{0\uparrow }={{\rm{\Omega }}}_{0{\rm{R}}}$. Despite these simplifying assumptions, the system is still very complex and gives rise to intriguing phenomena.

We find the stationary states of the system by self-consistently solving (4)–(10) in the mean-field regime in the parameter space {η, Δ}. This amounts to omitting quantum fluctuations and replacing the atomic and cavity field operators by their corresponding quantum averages: ${\hat{\psi }}_{\tau }\to {\psi }_{\tau }:= \langle {\hat{\psi }}_{\tau }\rangle$, ${\hat{a}}_{j}\to {\alpha }_{j}:= \langle {\hat{a}}_{j}\rangle$, and ${\hat{b}}_{j}\to {\beta }_{j}:= \langle {\hat{b}}_{j}\rangle$. The parameters η and Δ are related, respectively, to the intensity and the frequency of the external pump lasers and can be readily tuned in experiment.

The effective atomic Hamiltonian density (5) is λ/2 periodic. Nonetheless, solving the equations for different numbers of unit cells (of length λ/2) reveals that the atomic condensate wave functions ψ_τ(z) are λ periodic in parameter regimes possessing cavity-induced spin–orbit coupling. Therefore, we always solve the mean-field equations corresponding to (4)–(10) in two unit cells of total length λ with periodic boundary conditions. The relation between cavity-induced spin–orbit coupling and the doubling of the periodicity of the condensate wave functions will be discussed in more details in section 3.3.

3.1. Atomic phase diagram

The mean-field density-wave order parameters ${{ \mathcal N }}_{\tau }=\langle {\hat{{ \mathcal N }}}_{\tau }\rangle$, see (10a), can be used to characterize the density structure of each BEC component. They quantify the magnitude of the density modulations of each BEC component, where a zero density-wave order parameter corresponds to a homogeneous density distribution. Because of the symmetric choice of the parameters as described above, we always find that the absolute values of the two density-wave order parameters are equal to one another, $| {{ \mathcal N }}_{\downarrow }| =| {{ \mathcal N }}_{\uparrow }|$.

On the other hand, the mean-field local pseudospin vector ${\bf{s}}(z)=({s}_{x}(z),{s}_{y}(z),{s}_{z}(z))$ $=\left. \langle \left.{\psi }_{{\rm{eff}}}(z\right)\right|\left|\vec{\sigma }\right|\left|\left.{\psi }_{{\rm{eff}}}(z\right)\right. \rangle$, where $\left|{\psi }_{\mathrm{eff}}(z)\right\rangle := {\left({\psi }_{\uparrow }(z),{\psi }_{\downarrow }(z)\right)}^{T}$ and $\left|\vec{\sigma }\right|=({\sigma }_{x},{\sigma }_{y},{\sigma }_{z})$ is the vector of the Pauli matrices, can be used to illustrate the spatial spin texture of the steady-states. The z component of the local pseudospin ${s}_{z}(z)=[| {\psi }_{\uparrow }(z){| }^{2}-| {\psi }_{\downarrow }(z){| }^{2}]/2$ is zero everywhere in all parameter regimes due to the symmetric choice of the parameters. Therefore, the local pseudospin vector always lies in the x–y plane. We find that the transverse local pseudospin vector varies in space in all parameter regimes. In some regimes ${\bf{s}}(z)$ exhibits a λ/2-periodic 'spin-wave' of ferromagnetic-magnon nature, meaning that the spin angle $\phi (z):= \arctan ({s}_{y}(z)/{s}_{x}(z))$ only sweeps a small sector within the interval [0, π/2] over a λ/2 distance. While for other parameter regimes, ϕ(z) sweeps a full 2π angle in the x–y plane over a λ/2 distance, leading to a λ/2-periodic 'spin spiral' of topological Skyrmionic nature [62, 63]. As a result, the mean-field collective atomic pseudospins ${S}_{\pm }=\langle {\hat{S}}_{\pm }\rangle =\int [{s}_{x}(z)\,\pm {{\rm{i}}{s}}_{y}(z)]{\rm{d}}{z}$ and the spin-wave order parameters ${{ \mathcal S }}_{\pm }^{(\pm )}=\langle {\hat{{ \mathcal S }}}_{\pm }^{(\pm )}\rangle =\int {{\rm{e}}}^{\pm 2{\rm{i}}{k}{z}}[{s}_{x}(z)\pm {{\rm{i}}{s}}_{y}(z)]{\rm{d}}{z}$, see (10b) and (10c), exhibit different behaviors in the spin-wave and the spin-spiral states. The spin-wave order parameters ${{ \mathcal S }}_{-}^{(-)}={\left({{ \mathcal S }}_{+}^{+}\right)}^{\dagger }$ describe the dominant process where a photon is scattered from one pumped mode into the counterpropagating other pumped mode (${\hat{a}}_{+}\leftrightarrow {\hat{b}}_{-}$). In contrast, ${{ \mathcal S }}_{-}^{(+)}={\left({{ \mathcal S }}_{+}^{(-)}\right)}^{\dagger }$ are the order parameters for Raman scattering between the unpumped modes (${\hat{a}}_{-}\leftrightarrow {\hat{b}}_{+}$).

The spin-wave and spin-spiral states can be quantitatively distinguished by their distinct topological structures via an appropriate topological invariant. The relevant topological invariant to characterize the spin texture of the system is the winding number [65, 66]

$$\begin{eqnarray}&&{ \mathcal W }:= \displaystyle \frac{1}{2\pi }{\int }_{0}^{\lambda /2}[{\partial }_{z}\phi (z)]{\rm{d}}{z}=\displaystyle \frac{\phi (\lambda /2)-\phi (0)}{2\pi },\end{eqnarray} \tag{ 11 }$$

where ϕ(z) defines the direction of the local pseudospin vector s(z) in the x–y plane. Note that the angle ϕ(z) is tied to the relative phase between the two condensate wave functions ${\psi }_{\downarrow }(z)$ and ${\psi }_{\uparrow }(z)$. The winding number essentially counts the number of full rotations of the local pseudospin vector s(z) around the origin in one unit cell. Zero winding number corresponds to the topologically trivial spin-wave state, while a nonzero winding number (i.e. ${ \mathcal W }=+1$) indicates the topological spin-spiral state.

The atomic phase diagram of the system in the rescaled parameter space $\left\{\sqrt{N}\eta /{\omega }_{\mathrm{rec}},{\rm{\Delta }}/{\omega }_{\mathrm{rec}}\right\}$ is shown in figure 2 and displays three fundamentally different phases. In the first phase, which corresponds to the region below the solid red curve in the phase diagram, the density order parameters are nonzero, $| {{ \mathcal N }}_{\tau }| \ne 0$, while the winding number is zero, ${ \mathcal W }=0$. This implies that the atomic density distribution ${n}_{\tau }(z)=| {\psi }_{\tau }(z){| }^{2}$ has a (λ/2-periodic) crystalline order and the local pseudospin s(z) exhibits a (also λ/2-periodic) spin-wave texture. Hence, we refer to this phase as the DW-SW state. Recall that $| {{ \mathcal N }}_{\downarrow }| =| {{ \mathcal N }}_{\uparrow }|$ due to the symmetric choice of parameters. The amplitude of the density modulations increases with increasing η, indicated by the increasing density-wave order parameters $| {{ \mathcal N }}_{\tau }|$. Because of the direct pumping of two of the cavity modes and the presence of two-photon Raman processes which scatter photons into the unpumped modes without transfering momentum to the atoms as shown in figure 3(c), there is no threshold behavior for the onset of the density waves. This is in contrast to transversally pumped ring cavities [36, 45] and to the case where both pairs of modes $\left\{{\hat{a}}_{\pm },{\hat{b}}_{\pm }\right\}$ couple to the same atomic transition [59]. The angle in the x–y plane which is sweeped by the local pseudospin vector s(z) over a distance of λ/2, is always restricted to the interval [0, π/2] for the entire parameter regime of the DW-SW phase. Hence, the spin-wave oscillations remain small in this phase. A typical spin-wave texture in the DW-SW phase is illustrated in figure 4(a), where the local pseudospin vector s(z) exhibits small oscillations in the x–y plane. The change in the length of the local pseudospin vector is due to the density modulations of the condensates. Note that the local pseudospin vector s(z) always lies in the x–y plane, as s_z(z) = 0 owing to the symmetric choice of the parameters.

Zoom In Zoom Out Reset image size

The second phase, corresponding to the uniform black region in the phase diagram of figure 2, is the PW-SS state. In this regime the density-wave order parameters are identically zero, $| {{ \mathcal N }}_{\tau }| =0$, while the winding number is nonzero, ${ \mathcal W }=1$. Therefore, the condensate densities are homogeneous in this phase, while the local pseudospin exhibits a spin-spiral texture. Figure 4(b) depicts a representative spin-spiral state in this phase. As indicated by the winding number ${ \mathcal W }=1$, the local pseudospin vector s(z) sweeps a full 2π angle over a λ/2 distance. The length of the local pseudospin vector is constant in space due to the uniform condensate densities in this phase. Note that again the local pseudospin vector has no z component, s_z(z) = 0.

In the third phase, both the density order parameters and the winding number are nonzero, $\left\{| {{ \mathcal N }}_{\tau }| \ne 0,{ \mathcal W }=1\right\}$. This implies that the density wave and the spin-spiral coexist in this phase, hence the name DW-SS state. The BEC densities exhibit strong modulations in this phase. The local pseudospin vector ${\bf{s}}(z)$, therefore, inherits this and its length changes drastically in space, while again sweeping a full 2π angle over a λ/2 distance as shown in figure 4(c). As for the other two phases, the local pseudospin vector lies in the x–y plane in the DW-SS phase.

For small cavity detuning, ${\rm{\Delta }}\lesssim -9{\omega }_{\mathrm{rec}}$, the system becomes unstable (hashed region in figure 2). This is due to the fact that the effective relative cavity detuning, i.e. the dispersively shifted bare relative cavity detuning due to the presence of the atoms, becomes positive (i.e. blue detuned) resulting in cavity heating. In contrast to a single component BEC in a cavity where the dispersive shift is solely given by the factor U₀N, an analytical expression for the dispersive shift of the cavity detuning in our model is nontrivial due to the various dispersive terms and coupling terms in the effective cavity-field Hamiltonian density (8).

The collective atomic spins S_± and the spin-wave order parameters ${{ \mathcal S }}_{\pm }^{(\pm )}$ exhibit distinct behaviors in these three phases. Figures 5(a) and (b) show the absolute values of S_± (solid blue curves) and ${{ \mathcal S }}_{-}^{(-)}={\left({{ \mathcal S }}_{+}^{(+)}\right)}^{\dagger }$ (dashed red curves) as a function of the rescaled pump strength $\sqrt{N}\eta /{\omega }_{\mathrm{rec}}$ for constant cavity detunings ${\rm{\Delta }}/{\omega }_{\mathrm{rec}}=-10$ and −20, respectively. By increasing the pump strength η from zero, for ${\rm{\Delta }}=-10{\omega }_{\mathrm{rec}}$ in figure 5(a) the system undergoes a phase transition from the DW-SW state to the DW-SS state, while for ${\rm{\Delta }}=-20{\omega }_{\mathrm{rec}}$ in figure 5(b) the phase transition from the DW-SW to the DW-SS occurs indirectly via the intermediate PW-SS state (see figure 2). While the collective atomic spins S_± are nonzero in both DW-SW and DW-SS phases, it vanishes in the PW-SS phase. The latter can be understood by the fact that the local pseudospin vector s(z) has a constant length over space in the PW-SS phase and it does a full 2π rotation uniformly over one unit cell, resulting in ${S}_{\pm }=\int {s}_{x}(z){\rm{d}}{z}\pm \int {{\rm{i}}{s}}_{y}(z){\rm{d}}{z}=0$. On the other hand the spin-wave order parameters ${{ \mathcal S }}_{-}^{(-)}={\left({{ \mathcal S }}_{+}^{(+)}\right)}^{\dagger }$ are nonzero in all three phases, indicating spin modulations in all regimes.

Zoom In Zoom Out Reset image size

The winding number ${ \mathcal W }$ jumps from zero to one across the phase transitions from the DW-SW to the PW-SS/DW-SS, signaling that these are topological phase transitions (the red solid curve in the phase diagram in figure 2, which is somewhat ragged due to the extremely slow convergence of numerics around these phase boundaries). In addition, the density-wave order parameters ${{ \mathcal N }}_{\tau }$, the collective atomic spins S_±, and the spin-wave order parameters ${{ \mathcal S }}_{\pm }^{(\pm )}$ exhibit discontinuous behaviors on the onset of these phase transitions (see figures 2 and 5). This indicates that the topological phase transitions from the DW-SW state to the PW-SS/DW-SS states also have first-order characteristics.

Although the atomic parameters $\left\{{{ \mathcal N }}_{\tau },{S}_{\pm },{{ \mathcal S }}_{\pm }^{(\pm )}\right\}$ change continuously across the PW-SS to DW-SS phase transition (see figures 2 and 5(b)), they exhibit nonanalytical behavior, a characteristic of a second-order phase transition. Therefore, the phase transition from the PW-SS to the DW-SS is second-order (yellow dashed curve in the phase diagram in figure 2). Note that the winding number ${ \mathcal W }$ is one in both phases, and therefore it does not change across this phase transition.

The phase boundary between the DW-SW and PW-SS phases is linear, whereas the other phase boundaries show more complex behaviors. All phase boundaries with different natures (i.e. topological first-order, and topologically trivial second-order phase transitions) meet at a single tricritical point, denoted by a green dot in the phase diagram in figure 2.

3.2. Cavity-field phase diagram

The mean-field amplitudes of the unpumped cavity modes ${\alpha }_{-}$ and ${\beta }_{+}$ can be exploited as the cavity-field order parameters to further characterize the system. The inset of figure 2 shows the absolute values of the rescaled unpumped modes $| {\alpha }_{-}| /\sqrt{N}=| {\beta }_{+}| /\sqrt{N}$ in the rescaled parameter space $\left\{\sqrt{N}\eta /{\omega }_{\mathrm{rec}},{\rm{\Delta }}/{\omega }_{\mathrm{rec}}\right\}$. The absolute values of the unpumped modes are equal to each other once again due to the symmetric choice of the parameters. The cavity-field phase diagram has a similar form as the atomic phase diagram. In particular, the field phase boundaries coincide precisely with the atomic phase boundaries. The field order parameters $\left\{{\alpha }_{-},{\beta }_{+}\right\}$ are nonzero in the DW-SW and DW-SS phases and increase monotonically by increasing pump strength. However, they are identically zero in the entire PW-SS phase. This can be understood by noting that the density-order parameters ${{ \mathcal N }}_{\tau }$ and the collective spins S_± are zero in this regime as discussed above. The dynamics of the two umpumped modes $\left\{{\alpha }_{-},{\beta }_{+}\right\}$ then decouple completely from the pumped ones $\left\{{\alpha }_{+},{\beta }_{-}\right\}$ (see (8)), and, therefore, no photons are scattered into these unpumped modes in this phase. This leads to uniform potentials and Raman coupling (see (6) and (7)), which in turn results in homogeneous condensate densities in a self-consistent manner. This signifies the nonlinear dynamical nature of the system. At this point it should be mentioned that the mirror backscattering into the unpumped modes which is always observed in such types of setups [67] will lead to slight modifications of this phase diagram. Nevertheless, the qualitative picture will stay the same.

Figures 6(a) and (b) show cuts through the field phase diagram along the rescaled pump strength $\sqrt{N}\eta /{\omega }_{\mathrm{rec}}$ at constant cavity detunings ${\rm{\Delta }}/{\omega }_{\mathrm{rec}}=-10$ and −20, together with the corresponding winding numbers ${ \mathcal W }$. The field-order parameters $\left\{{\alpha }_{-},{\beta }_{+}\right\}$ exhibit similar behavior as the atomic parameters $\left\{{{ \mathcal N }}_{\tau },{S}_{\pm },{{ \mathcal S }}_{\pm }^{(\pm )}\right\}$. They display first-order (second-order) characteristics across the phase transitions from the DW-SW state to the PW-SS/DW-SS states (from the PW-SS state to the DW-SS phase), in accordance with the atomic phase transitions. Likewise, the first-oder phase transition coincides with the topological phase transition, where the winding ${ \mathcal W }$ jumps from zero to one. Therefore, there is a one to one correspondence between the atomic parameters $\left\{{{ \mathcal N }}_{\tau },{S}_{\pm },{{ \mathcal S }}_{\pm }^{(\pm )}\right\}$ on the one hand and the cavity-field order parameters $\left\{{\alpha }_{-},{\beta }_{+}\right\}$ on the other hand. As a consequence, all the quantum phase transitions (and their natures) can be mapped out nondestructively through the cavity outputs. This is an important and distinct feature of the system.

Zoom In Zoom Out Reset image size

3.3. Atomic momentum distributions and cavity-induced spin–orbit coupling

The discrete momentum exchange between the atoms and the light fields allows the decomposition of the condensate wave functions into plane waves ${\psi }_{\tau }(z)={\sum }_{j=-\infty }^{\infty }{c}_{\tau ,j}{{\rm{e}}}^{{\rm{i}}{j}{k}{z}}$. The absolute values of the probability amplitudes ${c}_{\tau ,j}$ of the lowest six momentum states $j\in \left\{0,\pm 1,\pm 2,+3\right\}$ in the rescaled parameter space $\left\{\sqrt{N}\eta /{\omega }_{\mathrm{rec}},{\rm{\Delta }}/{\omega }_{\mathrm{rec}}\right\}$ for each condensate component τ are shown in figure 7. Note that the even and odd momentum states do not coexist together. The boundary separating even and odd momenta coincides precisely with the topological phase boundary between the DW-SW state and the PW-SS/DW-SS states, illustrated in figure 2. The region where even (odd) momenta are occupied corresponds to the DW-SW (PW-SS or DW-SS) phase. In the DW-SW phase the zero momentum c_τ,0 is the dominant state for both condensates and the nonzero higher momenta result in density modulations. In the PW-SS phase, the condensate wave functions are solely composed of one momentum component $| {c}_{\downarrow ,1}| =| {c}_{\uparrow ,-1}| =1/\sqrt{2}$, as expected for a homogeneous condensate. While in the DW-SS phase, higher odd momenta are also populated, leading to density modulations. For the sake of clarity, vertical cuts along the rescaled pump strength $\sqrt{N}\eta /{\omega }_{\mathrm{rec}}$ for constant cavity detuning Δ = −20ω_rec of these momentum phase diagrams are also shown in figure 8.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The phase boundaries where the even momenta completely deplete and the odd momenta are populated coincide exactly with the first-order topological phase transitions from the DW-SW state to the PW-SS/DW-SS state. That is, the onset of the occupation of the odd momenta marks the appearance of spin spirals and the discrete jump of the winding number ${ \mathcal W }$ from zero to one; see figures 5(b), 6(b) and 8. This is intimately connected to the emergence of cavity-induced spin–orbit coupling for the atoms. The transition from the even to the odd momenta can be triggered and relaxed by (even very weak) two-body contact interactions, which are not included explicitly in our model.

The emergence of cavity-induced spin–orbit coupling can be most easily seen in the PW-SS phase, where the umpumped cavity modes $\left\{{\alpha }_{-},{\beta }_{+}\right\}$ are zero. The effective atomic Hamiltonian density (5) in the mean-field approximation then simplifies to

$$\begin{eqnarray}&&{{ \mathcal H }}_{\mathrm{SOC}}=\displaystyle \frac{1}{2m}{\left({{pI}}_{2\times 2}-{ \hbar }k{\sigma }_{z}\right)}^{2}+\displaystyle \frac{{ \hbar }{U}_{0}}{2}(| {\alpha }_{+}{| }^{2}-| {\beta }_{-}{| }^{2}){\sigma }_{z}+{ \hbar }{{\rm{\Omega }}}_{0{\rm{R}}}\left({\alpha }_{+}^{* }{\beta }_{-}{\sigma }_{\downarrow \uparrow }+{\alpha }_{+}{\beta }_{-}^{* }{\sigma }_{\uparrow \downarrow }\right),\end{eqnarray} \tag{ 12 }$$

after applying a unitary transformation [46]. Here ${I}_{2\times 2}$ is the identity matrix in the pseudospin space, σ_z is the third Pauli matrix and ${\sigma }_{\downarrow \uparrow }$ and ${\sigma }_{\uparrow \downarrow }$ are the transition matrices in the pseudospin basis. The Hamiltonian (12) has exactly the form of an equal Rashba–Dresselhaus spin–orbit coupled Hamiltonian, saving that the Raman coupling now depends on the cavity fields $\left\{{\alpha }_{+},{\beta }_{-}\right\}$ and is determined self-consistently. This Hamiltonian has been studied before in [46, 47, 49] and indeed exhibits characteristics of spin–orbit coupled quantum gases, with extra features resulting from the dynamical nature of the synthetic spin–orbit coupling.

The effect of the spin–orbit coupling can be seen in the momentum distributions of the condensate wavefuctions in the PW-SS phase, where different pseudospin states are coupled to different momentum states. That is, the pseudospin down is solely coupled to the $+{ \hbar }k$ momentum (recall that in the PW-SS phase $| {c}_{\downarrow ,1}| =1/\sqrt{2}$), while the the pseudospin up is only coupled to the $-{ \hbar }k$ momentum ($| {c}_{\uparrow ,-1}| =1/\sqrt{2}$). Since the PW-SS phase sets in at large pump strengths η, the effective Raman transition rate ${{\rm{\Omega }}}_{0{\rm{R}}}{\alpha }_{+}^{* }{\beta }_{-}$ is, therefore, always large. Hence, the single-particle energy dispersion of the spin–orbit coupled Hamiltonian (12) possesses a single minimum at p = 0, as expected for large Raman transition rates [12]. This is also the reason that both condensates have equal particle numbers, $| {c}_{\downarrow ,1}| =| {c}_{\uparrow ,-1}| =1/\sqrt{2}$, as the state at p = 0 has an equal contribution from the up and down components due to the symmetric choice of the parameters.

Despite of the fact that in the DW-SS phase the effective atomic Hamiltonian density (5) cannot be recast in the usual form of the equal Rashba–Dresselhaus spin–orbit-coupled Hamiltonian (12), the cavity-induced synthetic spin–orbit coupling still manifests itself in the momentum distributions of the condensate wave functions. Although higher odd momenta are also populated, the ±${ \hbar }$k momenta are still the dominant states, and different pseudospin states are strongly coupled to solely one of them, $| {c}_{\downarrow ,1}| =| {c}_{\uparrow ,-1}| \lesssim 1/\sqrt{2}$. This is in sharp contrast to the DW-SW phase, where both pseudospin states couple to the same zero momentum state, $| {c}_{\downarrow ,0}| =| {c}_{\uparrow ,0}| \lesssim 1/\sqrt{2}$, resulting in no cavity-induced spin–orbit coupling.

The period doubling of the condensate wave functions in the spin–orbit coupled PW-SS and DW-SS regimes, as it was mentioned at the beginning of section 3, can be understood through the momentum decomposition of the wave functions in these phases. In these spin–orbit coupled states, the condensate wave functions ${\psi }_{\tau }(z)={\sum }_{j=-\infty }^{\infty }{c}_{\tau ,2j+1}{{\rm{e}}}^{{\rm{i}}(2j+1){kz}}$ are composed of solely odd momenta and are, therefore, λ periodic. This is in contrast to the λ/2 periodicity of the condensate wave functions ${\psi }_{\tau }(z)={\sum }_{j=-\infty }^{\infty }{c}_{\tau ,2j}{{\rm{e}}}^{2{\rm{i}}{j}{k}{z}}$ in the DW-SW state (which are comprised of only even momenta) and the Hamiltonian density (5).

The unpumped cavity modes $\left\{{\hat{a}}_{-},{\hat{b}}_{+}\right\}$ play an important role in the emergence of cavity-induced spin–orbit coupling beyond a threshold for the pump-strength. This can be understood by re-examining the possible spin flipping processes in the system (see figure 3). The essential photon-scattering processes for the spin–orbit coupling are the ones depicted in figure 3(d), where photons are scattered between the modes ${\hat{a}}_{\pm }\leftrightarrow {\hat{b}}_{\mp }$ (via the atomic pseudospin flipping $\left.\left|\downarrow \right.,\rangle \right\rangle \leftrightarrow \left.\left|\uparrow \right.,\rangle \right\rangle$) and a ±2${ \hbar }$k momentum is transferred to the atom. Whereas the scattering processes shown in figure 3(c), where photons are scattered between the modes ${\hat{a}}_{\pm }\leftrightarrow {\hat{b}}_{\pm }$ (again via the atomic pseudospin flipping $\left.\left|\downarrow \right.,\rangle \right\rangle \leftrightarrow \left.\left|\uparrow \right.,\rangle \right\rangle$) without any momentum kick to the atom, are not vital for the spin–orbit coupling. However, triggering the former processes costs more energy than the latter ones due to the atomic kinetic energy gain. Therefore, for lower pump strengths (i.e. in the DW-SW phase) the spin flipping processes with no momentum kick to the atoms are energetically favored and are the dominant processes. The essential spin–orbit coupling processes become energetically favored and dominant beyond the pump-strength threshold on the onset of the PW-SS and DW-SS states, where the sum of the kinetic energies of the odd momenta becomes less than the corresponding even ones.

The interplay between spin flipping processes with and without momentum transfer to the atom can be seen by comparing the collective atomic spins S_± and the spin-wave order parameters ${{ \mathcal S }}_{\pm }^{(\pm )}$. Recall that S_± (${{ \mathcal S }}_{\pm }^{(\pm )}$) quantifies the spin flipping processes without (with ±2${ \hbar }$k) momentum kick to the atom. As can be seen from figure 5, the collective spins S_± are the dominant quantities in the DW-SW phase, while the spin-wave order parameters ${{ \mathcal S }}_{\pm }^{(\pm )}$ become dominant only in the PW-SS and DW-SS phases. Consequently, cavity-induced spin–orbit coupling emerges only in the PW-SS and DW-SS regimes. This is in sharp contrast to the free space spin–orbit-coupled BEC, where spin–orbit coupling emerges at an infinitesimal Raman frequency [12].

The features of the atomic momentum distribution shown in figure 8 serve as a direct experimental observation tool for the predicted phase transitions via time-of-flight measurements. In particular, the topologically nontrivial spin textures are directly related to the equal population of the $\pm { \hbar }k$ momenta for the two BEC components. This can be seen by writing down the momentum decomposition of the wave function of the PW-SS state given as ${\psi }_{\downarrow }(z)={c}_{\downarrow ,1}{{\rm{e}}}^{{\rm{i}}{k}{z}}$ and ${\psi }_{\uparrow }(z)={c}_{\uparrow ,-1}{{\rm{e}}}^{-{\rm{i}}{k}{z}}$. The relative, space dependent phase between those wave functions $\phi (z)=2{kz}$ defines the direction of the pseudospin vector at every position $z\in [0,\lambda /2]$. Therefore, the spin performs a full circle over a $\lambda /2$ distance leading to the spiral spin texture.

In addition, it should be mentioned that for the experimental observation of spin–orbit coupling in this system, the backscattering into the same mode due to imperfections of the cavity mirrors, has to be small. This can be achieved by either choosing an optimized cavity geometry or via an additional light field which destructively interferes with the backscattered light [68].

In order to check the stability of our mean-field results and to obtain a deeper understanding of the system, we calculate the collective excitations of the system above the mean-field steady-states. To this end, we linearize the Heisenberg equations of motion (4) for quantum fluctuations of both atomic condensate wave functions $\delta {\psi }_{\tau }(x,t)=\delta {\psi }_{\tau }^{(+)}(x){{\rm{e}}}^{-{\rm{i}}\omega t}+{\left[\delta {\psi }_{\tau }^{(-)}(x)\right]}^{* }{{\rm{e}}}^{{\rm{i}}{\omega }^{* }t}$, and field mode fluctuations $\delta {\alpha }_{\pm }(t)=\delta {\alpha }_{\pm }^{(+)}{{\rm{e}}}^{-{\rm{i}}\omega t}+{\left[\delta {\alpha }_{\pm }^{(-)}\right]}^{* }{{\rm{e}}}^{{\rm{i}}{\omega }^{* }t}$ and $\delta {\beta }_{\pm }(t)=\delta {\beta }_{\pm }^{(+)}{{\rm{e}}}^{-{\rm{i}}\omega t}+{\left[\delta {\beta }_{\pm }^{(-)}\right]}^{* }{{\rm{e}}}^{{\rm{i}}{\omega }^{* }t}$ around the mean-field stationary solutions $\left\{{\psi }_{0\tau }(x),{\alpha }_{0\pm },{\beta }_{0\pm }\right\}$. The linearized equations can be recast in matrix form,

$$\begin{eqnarray}&&\omega {\bf{f}}={{\bf{M}}}_{{\rm{B}}}{\bf{f}},\end{eqnarray} \tag{ 13 }$$

where ${\bf{f}}$ is a vector composed of the atomic condensate and the field-mode fluctuations $\left\{\delta {\psi }_{\tau }^{(\pm )},\delta {\alpha }_{\pm }^{(\pm )},\delta {\beta }_{\pm }^{(\pm )}\right\}$ and ${{\bf{M}}}_{{\rm{B}}}$ is a (nonhermitian) Bogoliubov matrix. We relegate the details to appendix B.

We numerically diagonalize the Bogoliubov matrix ${{\bf{M}}}_{{\rm{B}}}$ as a function of the pump strength to obtain the collective excitation spectrum $\omega (\eta )$. Figure 9 shows the real part of the five lowest lying positive-frequency excitations of the system as a function of the rescaled pump strength $\sqrt{N}\eta /{\omega }_{\mathrm{rec}}$ at a fixed cavity detuning ${\rm{\Delta }}=-20{\omega }_{\mathrm{rec}}$.

Zoom In Zoom Out Reset image size

In the DW-SW regime corresponding to $\sqrt{N}\eta \lesssim 27{\omega }_{\mathrm{rec}}$, there exists a gapless Goldstone mode, i.e. $\mathrm{Re}(\omega )=0$. This is associated with the spontaneously broken continuous screw-like symmetry of the system in this regime. The underlying symmetry breaking processes in our model are the fluctuations of the cavity fields due to the discreteness of cavity decay. This symmetry breaking is a consequence of the occupation of the even (in particular, the zero) momentum states in this phase, which leads to wave functions of the form ${\psi }_{\tau }(z)=({c}_{\tau ,0}+{c}_{\tau ,\pm 2}{{\rm{e}}}^{\pm 2{\rm{i}}{k}{z}}+...)$. A spatial translation of the wave functions ${\psi }_{\tau }(z)\to {\psi }_{\tau }(z+{\rm{\Delta }}z)$ cannot be compensated by the phase rotation of the condensate wave functions ${\psi }_{\downarrow }\to {\psi }_{\downarrow }{{\rm{e}}}^{-{\rm{i}}{k}{\rm{\Delta }}z}$ and ${\psi }_{\uparrow }\to {\psi }_{\uparrow }{{\rm{e}}}^{{\rm{i}}{k}{\rm{\Delta }}z}$, because of the occupation of the zero momentum state ${c}_{\tau ,0}\ne 0$. Hence, the condensate wave functions in the DW-SW regime are not invariant under the continuous symmetry group of the Hamiltonian and, therefore, they spontaneously break the continuous symmetry of the system. The appearance of this extra gapless Goldstone mode, in addition to the gapless phonon mode resulting from the spontaneous breaking of the internal gauge symmetry which is not shown here, indicates that the DW-SW state is a supersolid, i.e. a state with two spontaneously broken continuous symmetries [14, 29, 45]. A supersolid has the characteristics of both a crystal and a superfluid, that is, it is a state with a long-range periodic density order which can flow without experiencing any friction force.

For infinitely small values of η, the collective excitations at frequencies $\sim {\omega }_{\mathrm{rec}}$ are fourfold degenerate. These excitations correspond to condensate fluctuations of both BEC components at momenta $\pm { \hbar }k$. With increasing η, one of these branches approaches the zero energy. The point where this excitation branch touches zero coincides exactly with the mean-field critical pump strength for the phase transition from the DW-SW state to the PW-SS phase (dashed line in figure 9). As soon as the gap of the lowest $\pm { \hbar }k$ branch closes, the odd momenta become the lower energy states and the synthetic spin–orbit coupling emerges. Hence, it is the interplay between the two lowest even and odd excitation branches (solid red and dotted blue in figure 9) that leads to emergent spin–orbit coupling in this system.

In the PW-SS regime all excitation branches are doubly degenerate and gapped, indicating that the continuous symmetry of the system is not broken. This is intuitively obvious since the spin spiral with no density modulation perfectly respects the screw-like symmetry of the system. This can be seen most readily by the fact that the condensate wave functions ${\psi }_{\downarrow }(z)={c}_{\downarrow ,1}{{\rm{e}}}^{{\rm{i}}{k}{z}}$ and ${\psi }_{\uparrow }(z)={c}_{\uparrow ,-1}{{\rm{e}}}^{-{\rm{i}}{k}{z}}$ are invariant under the screw-like symmetry transformation: ${\psi }_{\downarrow }(z)\to {\psi }_{\downarrow }(z+{\rm{\Delta }}z){{\rm{e}}}^{-{\rm{i}}{k}{\rm{\Delta }}z}={\psi }_{\downarrow }(z)$ and ${\psi }_{\uparrow }(z)\to {\psi }_{\uparrow }(z+{\rm{\Delta }}z){{\rm{e}}}^{{\rm{i}}{k}{\rm{\Delta }}z}={\psi }_{\uparrow }(z)$.

With further increasing pump strength η the degenerate branches start to split up at the value of η which perfectly coincides with the mean-field critical pump strength for the phase transition from the PW-SS state to the DW-SS regime (horizontal dash dotted line). One of the excitation branches (dashed black curve) exhibits an exotic 'quasi' gapless-mode behavior in large pump strengths. This can be understood by examining the condensate wave functions in this regime. For the sake of simplicity we restrict our argument to one pseudospin, say, ${\psi }_{\downarrow }(z)=({c}_{\downarrow ,1}{{\rm{e}}}^{{\rm{i}}{k}{z}}+{c}_{\downarrow ,-1}{{\rm{e}}}^{-{\rm{i}}{k}{z}}+{c}_{\downarrow ,3}{{\rm{e}}}^{3{\rm{i}}{k}{z}}+...$). In general, this state is not invariant under the screw-like symmetry transformation ${\psi }_{\downarrow }(z)\to {\psi }_{\downarrow }(z+{\rm{\Delta }}z){{\rm{e}}}^{-{\rm{i}}{k}{\rm{\Delta }}z}\ne {\psi }_{\downarrow }(z)$. That said, at the onset of the DW-SS state ${c}_{\downarrow ,1}$ is the dominant probability amplitude in this expansion as can be seen from figure 7, i.e. $| {c}_{\downarrow ,1}| \gg \left\{| {c}_{\downarrow ,-1}| ,| {c}_{\downarrow ,3}| ,\,\cdots \right\}$. Therefore, this wave function can be approximated as ${\psi }_{\downarrow }(z)\simeq {c}_{\downarrow ,1}{{\rm{e}}}^{{\rm{i}}{k}{z}}$, which approximately preserves the screw-like symmetry of the system. For larger pump strengths, however, the higher momentum-state coefficients are not negligible anymore, resulting in a state which breaks the screw-like symmetry of the system. This leads to the appearance of a quasi gapless mode $\sim 0{\omega }_{\mathrm{rec}}$ at large η, which never reaches the zero energy exactly because ${c}_{\downarrow ,1}$ remains dominant throughout the entire regime. This is in stark contrast to the DW-SW state, where even in the onset of this phase at very small pump strengths the screw-like symmetry is broken completely. This is due to the dominant population of the zero momentum state $| {c}_{\tau ,0}| \gg \left\{| {c}_{\tau ,\pm 2}| ,| {c}_{\tau ,\pm 4}| ,\,\cdots \right\}$ in this entire regime (recall that the zero momentum state is the one that breaks the screw-like symmetry of the system in the DW-SW state).

It should be mentioned that imaginary parts of the collective excitations remain remarkably small for a dissipative system, i.e. $\kappa \ne 0$. This is a hint for the dynamical stability of the phases which is also confirmed by performing a real time evolution of the stationary state solutions.

We theoretically studied an effective two-component BEC inside a ring cavity, which possesses two pairs of nearly resonant running-wave modes with orthogonal polarizations. Our proposed model takes into account both atomic internal and external degrees of freedom, as well as the field amplitude and polarizations degrees of freedom. We predict that even in the simplest symmetric choice of parameters, the interplay between various degrees of freedom already results in novel phases and exotic quantum phase transitions of different natures. All the phases and the quantum phase transitions between them can be readily realized by solely tuning the frequencies and powers of the pump lasers, relevant parameters in cavity-QED experiments [27, 28, 64]. Remarkably, all the quantum phase transitions, including the topological one, can be monitored directly through the cavity outputs. The system, therefore, serves as an analog quantum simulator for studying topological phase transitions and dynamics of topologically nontrivial phases in a well controlled environment with in situ monitoring.

Our proposal can be implemented with minor modifications to state-of-the-art experiments in cavity QED [20, 21, 26, 27, 64, 67, 69, 70]. The main experimental challenges are the minimization of mirror backscattering, and inevitable birefringence at the mirrors. The former issue can be resolved by an additional light field which destructively interferes with the backscattered light fields, whereas birefringence may be reduced by choosing optimized mirror coatings. Additional physics may arise for asymmetric choices of the parameters as well as the inclusion of large two-body contact interactions. However, we leave the investigation of these interesting issues for future works.

We thank C Simon, F Piazza, and C Zimmermann for fruitful discussions. We acknowledge support by the Austrian Science Fund FWF through the Projects SFB FoQuS P13 and No. I1697-N27.

Here we demonstrate how the model given in (4) can be obtained from the single-particle Hamiltonian density (1). The single-particle Hamiltonian density (1) can be transferred into the rotating frame of the pump lasers through $\tilde{{ \mathcal H }}={ \mathcal U }{ \mathcal H }{{ \mathcal U }}^{\dagger }+{\rm{i}}{ \hbar }({\partial }_{t}{ \mathcal U }){{ \mathcal U }}^{\dagger }$ and exploiting the unitary transformation

$$\begin{eqnarray}&&{ \mathcal U }=\exp \left\{{\rm{i}}\left[({\hat{a}}_{+}^{\dagger }{\hat{a}}_{+}+{\hat{a}}_{-}^{\dagger }{\hat{a}}_{-}-{\hat{\sigma }}_{\downarrow \downarrow }){\tilde{\omega }}_{a}+({\hat{b}}_{+}^{\dagger }{\hat{b}}_{+}+{\hat{b}}_{-}^{\dagger }{\hat{b}}_{-}-{\hat{\sigma }}_{\uparrow \uparrow }){\tilde{\omega }}_{b}\right]t\right\}.\end{eqnarray} \tag{ A.1 }$$

The corresponding many-body Hamiltonian expressed in the formalism of second quantization then reads,

$$\begin{eqnarray}&&H=\int {\rm{d}}{z}{\hat{{\boldsymbol{\Psi }}}}^{\dagger }{ \mathcal M }\hat{{\boldsymbol{\Psi }}}-{ \hbar }\displaystyle \sum _{j=+,-}\left({{\rm{\Delta }}}_{a}{\hat{a}}_{j}^{\dagger }{\hat{a}}_{j}+{{\rm{\Delta }}}_{b}{\hat{b}}_{j}^{\dagger }{\hat{b}}_{j}\right)+i{ \hbar }\left[{\eta }_{+}{\hat{a}}_{+}^{\dagger }+{\eta }_{-}{\hat{b}}_{-}^{\dagger }-{\rm{h}}.{\rm{c}}.\right],\end{eqnarray} \tag{ A.2 }$$

where $\hat{{\boldsymbol{\Psi }}}(z)={\left({\hat{\psi }}_{\downarrow }(z),{\hat{\psi }}_{\uparrow }(z),{\hat{\psi }}_{e}(z)\right)}^{{\mathsf{T}}}$, and

$$\begin{eqnarray}{ \mathcal M }=\left(\begin{array}{ccc}\tfrac{{\hat{p}}^{2}}{2m}-\tfrac{{ \hbar }\delta }{2} & 0 & {{ \mathcal M }}_{13}\\ 0 & \tfrac{{\hat{p}}^{2}}{2m}+\tfrac{{ \hbar }\delta }{2} & {{ \mathcal M }}_{23}\\ {{ \mathcal M }}_{31} & {{ \mathcal M }}_{32} & \tfrac{{\hat{p}}^{2}}{2m}-\tfrac{{ \hbar }}{2}({{\rm{\Delta }}}_{\downarrow }+{{\rm{\Delta }}}_{\uparrow })\end{array}\right),\end{eqnarray} \tag{ A.3 }$$

with the elements

$$\begin{eqnarray}\begin{array}{rcl}{{ \mathcal M }}_{31} & = & {{ \mathcal M }}_{13}^{\dagger }={ \hbar }{{ \mathcal G }}_{\downarrow }\left({{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{a}}_{+}+{{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{a}}_{-}\right),\\ {{ \mathcal M }}_{32} & = & {{ \mathcal M }}_{23}^{\dagger }={ \hbar }{{ \mathcal G }}_{\uparrow }\left({{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{b}}_{+}+{{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{b}}_{-}\right).\end{array}\end{eqnarray} \tag{ A.4 }$$

The constant term $[({\omega }_{\downarrow }+{\tilde{\omega }}_{a})/2+({\omega }_{\uparrow }+{\tilde{\omega }}_{b})/2]{I}_{3\times 3}$ is omitted. The dynamics of the atomic and cavity field operators can be determined by simultaneously solving the following Heisenberg equations of motion ${\rm{i}}{ \hbar }{\partial }_{t}{\hat{\psi }}_{\tau }=[{\psi }_{\tau },H]$ and ${\rm{i}}{ \hbar }{\partial }_{t}{\hat{a}}_{j}/{\hat{b}}_{j}=[{\hat{a}}_{j}/{\hat{b}}_{j},H]-{\rm{i}}{ \hbar }\kappa {\hat{a}}_{j}/{\hat{b}}_{j}$. Substituting the Hamiltonian (A.2) leads to the following set of coupled differential equations

$$\begin{eqnarray}\begin{array}{rcl}{\rm{i}}{ \hbar }\tfrac{\partial }{\partial t}{\hat{\psi }}_{\downarrow } & = & \left(\tfrac{{\hat{p}}^{2}}{2m}-\tfrac{{ \hbar }}{2}\delta \right){\hat{\psi }}_{\downarrow }+{ \hbar }{{ \mathcal G }}_{\downarrow }^{* }\left({{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{a}}_{+}^{\dagger }+{{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{a}}_{-}^{\dagger }\right){\hat{\psi }}_{{e}},\\ {\rm{i}}{ \hbar }\tfrac{\partial }{\partial t}{\hat{\psi }}_{\uparrow } & = & \left(\tfrac{{\hat{p}}^{2}}{2m}+\tfrac{{ \hbar }}{2}\delta \right){\hat{\psi }}_{\uparrow }+{ \hbar }{{ \mathcal G }}_{\uparrow }^{* }\left({{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{b}}_{+}^{\dagger }+{{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{b}}_{-}^{\dagger }\right){\hat{\psi }}_{{e}},\\ {\rm{i}}{ \hbar }\tfrac{\partial }{\partial t}{\hat{a}}_{+} & = & -{ \hbar }({{\rm{\Delta }}}_{a}+{\rm{i}}\kappa ){\hat{a}}_{+}+{ \hbar }{{ \mathcal G }}_{\downarrow }^{* }\displaystyle \int {{\rm{d}}{z}{e}}^{-{\rm{i}}{k}{z}}{\hat{\psi }}_{\downarrow }^{\dagger }{\hat{\psi }}_{{e}}+{\rm{i}}{ \hbar }{\eta }_{+},\\ {\rm{i}}{ \hbar }\tfrac{\partial }{\partial t}{\hat{a}}_{-} & = & -{ \hbar }({{\rm{\Delta }}}_{a}+{\rm{i}}\kappa ){\hat{a}}_{-}+{ \hbar }{{ \mathcal G }}_{\downarrow }^{* }\displaystyle \int {{\rm{d}}{z}{e}}^{{\rm{i}}{k}{z}}{\hat{\psi }}_{\downarrow }^{\dagger }{\hat{\psi }}_{{e}},\\ {\rm{i}}{ \hbar }\tfrac{\partial }{\partial t}{\hat{b}}_{+} & = & -{ \hbar }({{\rm{\Delta }}}_{b}+{\rm{i}}\kappa ){\hat{b}}_{+}+{ \hbar }{{ \mathcal G }}_{\uparrow }^{* }\displaystyle \int {{\rm{d}}{z}{e}}^{-{\rm{i}}{k}{z}}{\hat{\psi }}_{\uparrow }^{\dagger }{\hat{\psi }}_{{e}},\\ {\rm{i}}{ \hbar }\tfrac{\partial }{\partial t}{\hat{b}}_{-} & = & -{ \hbar }({{\rm{\Delta }}}_{b}+{\rm{i}}\kappa ){\hat{b}}_{-}+{ \hbar }{{ \mathcal G }}_{\uparrow }^{* }\displaystyle \int {{\rm{d}}{z}{e}}^{{\rm{i}}{k}{z}}{\hat{\psi }}_{\uparrow }^{\dagger }{\hat{\psi }}_{{e}}+{\rm{i}}{ \hbar }{\eta }_{-}.\end{array}\end{eqnarray} \tag{ A.5 }$$

If the relative atomic detunings ${{\rm{\Delta }}}_{\uparrow }$ and ${{\rm{\Delta }}}_{\downarrow }$ are large compared to the two-photon detuning δ and the atom-photon couplings ${{ \mathcal G }}_{\uparrow }$ and ${{ \mathcal G }}_{\downarrow }$, the atomic excited state reaches a steady-state on a short time scale and its dynamics can be eliminated adiabatically. By setting ${\partial }_{t}{\hat{\psi }}_{e}=0$ in the Heisenberg equation of motion for ${\hat{\psi }}_{e}$ and omitting the kinetic energy compared to the term proportional to ${ \hbar }({{\rm{\Delta }}}_{\uparrow }+{{\rm{\Delta }}}_{\downarrow })$, we obtain the steady-state field operator of the atomic excited state

$$\begin{eqnarray}&&{\hat{\psi }}_{{e}}^{\mathrm{ss}}\simeq \displaystyle \frac{2}{({{\rm{\Delta }}}_{\downarrow }+{{\rm{\Delta }}}_{\uparrow })}\left[{{ \mathcal G }}_{\downarrow }\left({{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{a}}_{+}+{{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{a}}_{-}\right){\hat{\psi }}_{\downarrow }+{{ \mathcal G }}_{\uparrow }\left({{\rm{e}}}^{{\rm{i}}{k}{z}}{\hat{b}}_{+}+{{\rm{e}}}^{-{\rm{i}}{k}{z}}{\hat{b}}_{-}\right){\hat{\psi }}_{\uparrow }\right].\end{eqnarray} \tag{ A.6 }$$

This steady-state field operator can be substituted into the Heisenberg equations of motion (A.5) resulting in a set of six coupled nonlinear equations for $\left\{{\hat{\psi }}_{\downarrow },{\hat{\psi }}_{\uparrow },{\hat{a}}_{\pm },{\hat{b}}_{\pm }\right\}$ given in (4).

In this appendix we describe the calculations leading to the collective excitation spectrum presented in section 4 in more detail. Plugging the ansatz ${\psi }_{\tau }(x,t)={\psi }_{0\tau }(x)+\delta \psi (x,t)$, ${\alpha }_{\pm }(t)={\alpha }_{\pm 0}+\delta {\alpha }_{\pm }(t)$ and ${\beta }_{\pm }(t)={\beta }_{\pm 0}+\delta {\beta }_{\pm }(t)$ into the mean field version of (4) and performing the Bogoliubov transformation as it is already discussed in section 4 leads to the following linearized equations for the modes

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\alpha }_{+}^{(+)} & = & -{\tilde{{\rm{\Delta }}}}_{a}\delta {\alpha }_{+}^{(+)}+{U}_{0\downarrow }{{ \mathcal N }}_{\downarrow }\delta {\alpha }_{+}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{S}_{-}\delta {\beta }_{+}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal S }}_{-}^{(1/2)}\delta {\beta }_{-}^{(+)}\\ & & +{U}_{0\downarrow }{{ \mathcal A }}_{+* }^{\downarrow }\delta {\psi }_{\downarrow }^{(+)}+{U}_{0\downarrow }{{ \mathcal A }}_{+}^{\downarrow }\delta {\psi }_{\downarrow }^{(-)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal B }}_{+* }^{\downarrow }\delta {\psi }_{\uparrow }^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal B }}_{+}^{\uparrow }\delta {\psi }_{\downarrow }^{(-)},\end{array}\end{eqnarray} \tag{ B.1 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\alpha }_{+}^{(-)} & = & {\tilde{{\rm{\Delta }}}}_{a}^{* }\delta {\alpha }_{+}^{(-)}-{U}_{0\downarrow }{{ \mathcal N }}_{\downarrow }^{* }\delta {\alpha }_{-}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{S}_{-}^{* }\delta {\beta }_{+}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal S }}_{-}^{(1)* }\delta {\beta }_{-}^{(-)}\\ & & -{U}_{0\downarrow }{{ \mathcal A }}_{+* }^{\downarrow * }\delta {\psi }_{\downarrow }^{(-)}-{U}_{0\downarrow }{{ \mathcal A }}_{+}^{\downarrow * }\delta {\psi }_{\downarrow }^{(+)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal B }}_{+* }^{\downarrow * }\delta {\psi }_{\uparrow }^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal B }}_{+}^{\uparrow * }\delta {\psi }_{\downarrow }^{(+)},\end{array}\end{eqnarray} \tag{ B.2 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\alpha }_{-}^{(+)} & = & -{\tilde{{\rm{\Delta }}}}_{a}\delta {\alpha }_{-}^{(+)}+{U}_{0\downarrow }{{ \mathcal N }}_{\downarrow }^{* }\delta {\alpha }_{+}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{S}_{-}\delta {\beta }_{-}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal S }}_{-}^{(2)}\delta {\beta }_{+}^{(+)}\\ & & +{U}_{0\downarrow }{{ \mathcal A }}_{-* }^{\downarrow }\delta {\psi }_{\downarrow }^{(+)}+{U}_{0\downarrow }{{ \mathcal A }}_{-}^{\downarrow }\delta {\psi }_{\downarrow }^{(-)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal B }}_{-* }^{\downarrow }\delta {\psi }_{\uparrow }^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal B }}_{-}^{\uparrow }\delta {\psi }_{\downarrow }^{(-)},\end{array}\end{eqnarray} \tag{ B.3 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\alpha }_{-}^{(-)} & = & {\tilde{{\rm{\Delta }}}}_{a}^{* }\delta {\alpha }_{-}^{(-)}-{U}_{0\downarrow }{{ \mathcal N }}_{\downarrow }\delta {\alpha }_{+}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{S}_{-}^{* }\delta {\beta }_{-}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal S }}_{-}^{(2)* }\delta {\beta }_{+}^{(-)}\\ & & -{U}_{0\downarrow }{{ \mathcal A }}_{-* }^{\downarrow * }\delta {\psi }_{\downarrow }^{(-)}-{U}_{0\downarrow }{{ \mathcal A }}_{-}^{\downarrow * }\delta {\psi }_{\downarrow }^{(+)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal B }}_{-* }^{\downarrow * }\delta {\psi }_{\uparrow }^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal B }}_{-}^{\uparrow * }\delta {\psi }_{\downarrow }^{(+)},\end{array}\end{eqnarray} \tag{ B.4 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\beta }_{+}^{(+)} & = & -{\tilde{{\rm{\Delta }}}}_{a}\delta {\beta }_{+}^{(+)}+{U}_{0\uparrow }{{ \mathcal N }}_{\uparrow }\delta {\beta }_{-}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{S}_{-}^{* }\delta {\alpha }_{+}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal S }}_{-}^{(2)* }\delta {\alpha }_{-}^{(+)}\\ & & +{U}_{0\uparrow }{{ \mathcal B }}_{+* }^{\uparrow }\delta {\psi }_{\uparrow }^{(+)}+{U}_{0\uparrow }{{ \mathcal B }}_{+}^{\uparrow }\delta {\psi }_{\uparrow }^{(-)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal A }}_{+* }^{\uparrow }\delta {\psi }_{\downarrow }^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal A }}_{+}^{\downarrow }\delta {\psi }_{\uparrow }^{(-)},\end{array}\end{eqnarray} \tag{ B.5 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\beta }_{+}^{(-)} & = & {\tilde{{\rm{\Delta }}}}_{a}^{* }\delta {\beta }_{+}^{(-)}-{U}_{0\uparrow }{{ \mathcal N }}_{\uparrow }^{* }\delta {\beta }_{-}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}{S}_{-}\delta {\alpha }_{+}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal S }}_{-}^{(2)}\delta {\alpha }_{-}^{(-)}\\ & & -{U}_{0\uparrow }{{ \mathcal B }}_{+* }^{\uparrow * }\delta {\psi }_{\uparrow }^{(-)}-{U}_{0\uparrow }{{ \mathcal B }}_{+}^{\uparrow * }\delta {\psi }_{\uparrow }^{(+)}-{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal A }}_{+* }^{\uparrow * }\delta {\psi }_{\downarrow }^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal A }}_{+}^{\downarrow * }\delta {\psi }_{\uparrow }^{(+)},\end{array}\end{eqnarray} \tag{ B.6 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\beta }_{-}^{(+)} & = & -{\tilde{{\rm{\Delta }}}}_{a}\delta {\beta }_{-}^{(+)}+{U}_{0\uparrow }{{ \mathcal N }}_{\uparrow }^{* }\delta {\beta }_{+}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{S}_{-}^{* }\delta {\alpha }_{-}^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal S }}_{-}^{(1)* }\delta {\alpha }_{+}^{(+)}\\ & & +{U}_{0\uparrow }{{ \mathcal B }}_{-* }^{\uparrow }\delta {\psi }_{\uparrow }^{(+)}+{U}_{0\uparrow }{{ \mathcal B }}_{-}^{\uparrow }\delta {\psi }_{\uparrow }^{(-)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal A }}_{-* }^{\uparrow }\delta {\psi }_{\downarrow }^{(+)}+{{\rm{\Omega }}}_{0{\rm{R}}}^{* }{{ \mathcal A }}_{-}^{\downarrow }\delta {\psi }_{\uparrow }^{(-)},\end{array}\end{eqnarray} \tag{ B.7 }$$

$$\begin{eqnarray}&&\omega \delta {\beta }_{-}^{(-)}={\tilde{{\rm{\Delta }}}}_{a}^{* }\delta {\beta }_{-}^{(-)}-{U}_{0\uparrow }{{ \mathcal N }}_{\uparrow }\delta {\beta }_{+}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}{S}_{-}\delta {\alpha }_{-}^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal S }}_{-}^{(1)}\delta {\alpha }_{+}^{(-)}-{U}_{0\uparrow }{{ \mathcal B }}_{-* }^{\uparrow * }\delta {\psi }_{\uparrow }^{(-)}\end{eqnarray} \tag{ B.8 }$$

$$\begin{eqnarray}&&-{U}_{0\uparrow }{{ \mathcal B }}_{-}^{\uparrow * }\delta {\psi }_{\uparrow }^{(+)}-{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal A }}_{-* }^{\uparrow * }\delta {\psi }_{\downarrow }^{(-)}-{{\rm{\Omega }}}_{0{\rm{R}}}{{ \mathcal A }}_{-}^{\downarrow * }\delta {\psi }_{\uparrow }^{(+)},\end{eqnarray} \tag{ B.9 }$$

where we introduced the following shorthand notations

$$\begin{eqnarray}\begin{array}{rcl}{{ \mathcal A }}_{\pm }^{\downarrow \uparrow }\xi : & = & \displaystyle \int {A}_{\pm }^{\downarrow \uparrow }\xi {\rm{d}}{z},\\ {{ \mathcal A }}_{\pm * }^{\downarrow \uparrow }\xi : & = & \displaystyle \int {A}_{\pm * }^{\downarrow \uparrow }\xi {\rm{d}}{z},\\ {{ \mathcal B }}_{\pm }^{\downarrow \uparrow }\xi : & = & \displaystyle \int {B}_{\pm }^{\downarrow \uparrow }\xi {\rm{d}}{z},\\ {{ \mathcal B }}_{\pm * }^{\downarrow \uparrow }\xi : & = & \displaystyle \int {B}_{\pm * }^{\downarrow \uparrow }\xi {\rm{d}}{z},\end{array}\end{eqnarray} \tag{ B.10 }$$

with

$$\begin{eqnarray}\begin{array}{rcl}{A}_{\pm }^{\downarrow \uparrow } & = & {\psi }_{0}^{\downarrow \uparrow }\left({\alpha }_{0\pm }+{{\rm{e}}}^{\mp 2{\rm{i}}{k}{z}}{\alpha }_{0\mp }\right),\\ {A}_{\pm * }^{\downarrow \uparrow } & = & {\psi }_{0}^{\downarrow \uparrow * }\left({\alpha }_{0\pm }+{{\rm{e}}}^{\mp 2{\rm{i}}{k}{z}}{\alpha }_{0\mp }\right),\\ {B}_{\pm }^{\downarrow \uparrow } & = & {\psi }_{0}^{\downarrow \uparrow }\left({\beta }_{0\pm }+{{\rm{e}}}^{\mp 2{\rm{i}}{k}{z}}{\beta }_{0\mp }\right),\\ {B}_{\pm * }^{\downarrow \uparrow } & = & {\psi }_{0}^{\downarrow \uparrow * }\left({\beta }_{0\pm }+{{\rm{e}}}^{\mp 2{\rm{i}}{k}{z}}{\beta }_{0\mp }\right).\end{array}\end{eqnarray} \tag{ B.11 }$$

The linearized equations for the atomic degrees of freedom read

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\psi }_{\downarrow }^{(+)} & = & \tfrac{1}{{ \hbar }}\left[{D}_{\downarrow ,1}-\mu \right]\delta {\psi }_{\downarrow }^{(+)}+{ \hbar }{{\rm{\Omega }}}_{{\rm{R}}}(z)\delta {\psi }_{\uparrow }^{(+)}\\ & & +{U}_{0\downarrow }\left({A}_{+* }^{\downarrow * }\delta {\alpha }_{+}^{(+)}+{A}_{-* }^{\downarrow * }\delta {\alpha }_{-}^{(+)}+{A}_{+}^{\downarrow }\delta {\alpha }_{+}^{(-)}+{A}_{-}^{\downarrow }\delta {\alpha }_{-}^{(-)}\right)\\ & & +{{\rm{\Omega }}}_{0{\rm{R}}}\left({B}_{+}^{\uparrow }\delta {\alpha }_{+}^{(-)}+{B}_{-}^{\uparrow }\delta {\alpha }_{-}^{(-)}+{A}_{+* }^{\uparrow * }\delta {\beta }_{+}^{(+)}+{A}_{-* }^{\uparrow * }\delta {\beta }_{-}^{(+)}\right),\end{array}\end{eqnarray} \tag{ B.12 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\psi }_{\downarrow }^{(-)} & = & -\tfrac{1}{{ \hbar }}{\left[{D}_{\downarrow ,1}-\mu \right]}^{* }\delta {\psi }_{\downarrow }^{(-)}-{ \hbar }{{\rm{\Omega }}}_{{\rm{R}}}^{* }(z)\delta {\psi }_{\uparrow }^{(-)}\\ & & -{U}_{0\downarrow }\left({A}_{+* }^{\downarrow }\delta {\alpha }_{+}^{(-)}+{A}_{-* }^{\downarrow }\delta {\alpha }_{-}^{(-)}+{A}_{+}^{\downarrow * }\delta {\alpha }_{+}^{(+)}+{A}_{-}^{\downarrow * }\delta {\alpha }_{-}^{(+)}\right)\\ & & -{{\rm{\Omega }}}_{0{\rm{R}}}^{* }\left({B}_{+}^{\uparrow * }\delta {\alpha }_{+}^{(+)}+{B}_{-}^{\uparrow * }\delta {\alpha }_{-}^{(+)}+{A}_{+* }^{\uparrow }\delta {\beta }_{+}^{(-)}+{A}_{-* }^{\uparrow }\delta {\beta }_{-}^{(-)}\right),\end{array}\end{eqnarray} \tag{ B.13 }$$

$$\begin{eqnarray}\begin{array}{rcl}\omega \delta {\psi }_{\uparrow }^{(+)} & = & \tfrac{1}{{ \hbar }}\left[{D}_{\uparrow ,2}-\mu \right]\delta {\psi }_{\uparrow }^{(+)}+{ \hbar }{{\rm{\Omega }}}_{{\rm{R}}}^{* }(z)\delta {\psi }_{\downarrow }^{(+)}\\ & & +{U}_{0\uparrow }\left({B}_{+* }^{\uparrow * }\delta {\beta }_{+}^{(+)}+{B}_{-* }^{\uparrow * }\delta {\beta }_{-}^{(+)}+{B}_{+}^{\uparrow }\delta {\beta }_{+}^{(-)}+{B}_{-}^{\uparrow }\delta {\beta }_{-}^{(-)}\right)\\ & & +{{\rm{\Omega }}}_{0{\rm{R}}}^{* }\left({B}_{+* }^{\downarrow * }\delta {\alpha }_{+}^{(+)}+{B}_{-* }^{\downarrow * }\delta {\alpha }_{-}^{(+)}+{A}_{+}^{\downarrow }\delta {\beta }_{+}^{(-)}+{A}_{-}^{\downarrow }\delta {\beta }_{-}^{(-)}\right),\end{array}\end{eqnarray} \tag{ B.14 }$$

$$\begin{eqnarray}&&\omega \delta {\psi }_{\uparrow }^{(-)}=-\tfrac{1}{{ \hbar }}{\left[{D}_{\uparrow ,2}-\mu \right]}^{* }\delta {\psi }_{\uparrow }^{(-)}-{ \hbar }{{\rm{\Omega }}}_{{\rm{R}}}(z)\delta {\psi }_{\downarrow }^{(-)}\end{eqnarray} \tag{ B.15 }$$

$$\begin{eqnarray}\begin{array}{rcl} & & -{U}_{0\uparrow }\left({B}_{+* }^{\uparrow }\delta {\beta }_{+}^{(-)}+{B}_{-* }^{\uparrow }\delta {\beta }_{-}^{(-)}+{B}_{+}^{\uparrow * }\delta {\beta }_{+}^{(+)}+{B}_{-}^{\uparrow * }\delta {\beta }_{-}^{(+)}\right)\\ & & -{{\rm{\Omega }}}_{0{\rm{R}}}\left({B}_{+* }^{\downarrow }\delta {\alpha }_{+}^{(-)}+{B}_{-* }^{\downarrow }\delta {\alpha }_{-}^{(-)}+{A}_{+}^{\downarrow * }\delta {\beta }_{+}^{(+)}+{A}_{-}^{\downarrow * }\delta {\beta }_{-}^{(+)}\right),\end{array}\end{eqnarray} \tag{ B.16 }$$

where we introduce the shorthand notation ${D}_{\downarrow \uparrow ,i}:= -\displaystyle \frac{{p}^{2}}{2m}+{ \hbar }{U}_{\downarrow \uparrow }(z)+{\left(-1\right)}^{i}\tfrac{{ \hbar }\delta }{2}$ and μ denotes the chemical potential.

The set of equations (B.1) and (B.12) can be written in matrix form which results in (13). Since we do not have an analytical steady-state solution for the condensate wave functions and the cavity modes, we numerically diagonalize the Bogoliubov matrix for a numerically determined steady-state solution to obtain the collective excitation spectrum presented in figure 9.