Skip to main content

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

Cavity-induced quantum spin liquids


Quantum spin liquids provide paradigmatic examples of highly entangled quantum states of matter. Frustration is the key mechanism to favor spin liquids over more conventional magnetically ordered states. Here we propose to engineer frustration by exploiting the coupling of quantum magnets to the quantized light of an optical cavity. The interplay between the quantum fluctuations of the electro-magnetic field and the strongly correlated electrons results in a tunable long-range interaction between localized spins. This cavity-induced frustration robustly stabilizes spin liquid states, which occupy an extensive region in the phase diagram spanned by the range and strength of the tailored interaction. This occurs even in originally unfrustrated systems, as we showcase for the Heisenberg model on the square lattice.


Quantum spin liquids (QSLs) represent strongly correlated phases of matter, which are characterized by quantum fluctuations so dominant as to suppress magnetic ordering down to the lowest temperatures. Yet, the spins may be quantum mechanically entangled over long distances1,2,3. In Nature, QSLs are expected to occur in proximity to magnetic phases, but their existence often remains elusive. The key ingredient behind quantum spin liquid formation is, however, clearly identified: it is the presence of strong frustration, which disallows magnetic symmetry breaking, but need not be averse to, e.g., quantum mechanical singlet ordering. The routes towards frustration are manifold: one promising avenue is the focus on materials where magnetic ordering is penalized by the geometry of the lattice, such as for triangular, Kagomé or pyrochlore lattices4,5,6,7. Another one proceeds via the energetic competition of couplings of different ranges, like in the antiferromagnetic (AFM) J1J2 Heisenberg model or dipolar-interacting systems8,9,10, where the simultaneous appearance of nearest- and beyond-nearest-neighbour couplings counteracts global antiferromagnetism.

The challenge is then out to engineer robust QSL states of quantum condensed matter. Here, we will achieve this task by coupling an ordinary Heisenberg antiferromagnet on a square lattice to the electromagnetic field of an optical cavity.

The physical mechanism stabilizing the QSL takes the second route towards strong frustration to the extreme, by considering long-range AFM interactions described by an algebraically decaying spin–spin interaction ~rα including the case of all-to-all couplings α = 0, mediated by the cavity, cf. Fig. 1a. For the limiting case α = 0 and a cavity-induced interaction γ dominating over the nearest-neighbour Heisenberg coupling J, J/γ = 0, this realizes a state with long-range correlations mediated by singlets of arbitrarily large size (LRS). Away from this limit, and for decay exponents α 1, within a Schwinger–Boson approach, we find that the frustration imprinted by the cavity creates an extensive regime of QSL states. It is characterized by the absence of spontaneous symmetry breaking, and fractional excitations of both of a gapped (SL-I) and of gapless (SL-II) nature, cf. Fig. 1b. As a consequence of the underlying long-ranged interactions, correlations decay algebraically in both these phases.

Fig. 1: Implementation of a cavity-induced quantum spin liquid and phase diagram.

a Setup: a two-dimensional material, with nearest-neighbour exchange interaction J, is coupled to a cavity with fundamental frequencies ω,, whose field is represented by the light-blue arrows. The system is driven by an external laser with frequency ωL. b Level scheme: the electronic orbitals \(|{b}_{1,2}\rangle\), with energies ~ϵ1,2 are coupled to the auxiliary band \(\left|{b}_{3}\right\rangle\), with energy ~ϵ3 via the laser with Rabi frequency ΩL and the cavity modes a,. The third band is detuned from the laser by Δ3, and from the cavity modes by Δ,. c Phase diagram for the ground state of Hamiltonian (1), obtained from the bosonic spinon decomposition, as a function of the exponent α and of the coupling ratio γ/J, featuring spin-liquid (SL), long-range-single (LRS) and antiferromagnetic (AFM) phases. Error bars on the phase boundaries are within the symbols’ size. The inset shows the square lattice and the reciprocal one, with the respective primitive vectors.

In terms of physical implementation, we draw motivation from recent developments exploring the interplay of quantum materials with quantized light. This idea has been researched in the context of weakly correlated systems, mainly as a tool to reinforce superconductivity and other coherent many-body phases11,12,13,14,15,16,17,18,19,20. First works have also addressed the strong coupling regime, showing how existing phases can be manipulated in this way12,21,22. Here, we demonstrate that the coupling to a cavity can even induce phases that are not present in its absence: an unfrustrated AFM system is turned into a quantum spin liquid, provided the AFM interaction mediated by the cavity is sufficiently long-ranged and strong. To achieve these requirements, we develop a solid-state implementation harnessing localized electronic orbitals as effective spin degrees of freedom, coupled to the cavity modes via additional coherent laser drive, cf. Fig. 1a. This gives rise to quantum mechanically fluctuating, effective magnetic fields in all linearly independent spatial directions, which vanish on average. They thus counteract dynamically magnetization in any direction, but do not suppress the spin-singlet ordering, crucial for QSL states.



We consider a long-range SU(2)-symmetric Heisenberg model on a square lattice

$$H=J\mathop{\sum}\limits_{\langle i,j\rangle }\ {{{{{{{{\bf{S}}}}}}}}}_{i}\cdot {{{{{{{{\bf{S}}}}}}}}}_{j}+\gamma \mathop{\sum}\limits_{i\ne j}\ \frac{{{{{{{{{\bf{S}}}}}}}}}_{i}\cdot {{{{{{{{\bf{S}}}}}}}}}_{j}}{| {{{{{{{{\bf{r}}}}}}}}}_{ij}{| }^{\alpha }},$$

with \({{{{{{{{\bf{S}}}}}}}}}_{i}=({S}_{i}^{x},{S}_{i}^{y},{S}_{i}^{z})\) spin-1/2 operators on the lattice site i, J > 0 the nearest-neighbour AFM exchange, γ > 0 the strength of the long-range interaction modulated by the exponent α and rij ≡ ri − rj. Periodic boundary conditions are assumed. Before analysing the ground-state phase diagram of the Hamiltonian (1), let us qualitatively discuss the expected phases, starting with some known limiting cases. For γ = 0, the ground state of the Hamiltonian (1) displays Néel-like order23. For α = 0 and γJ, the long-range Hamiltonian is proportional to the total spin \({({\sum }_{i}{{{{{{{{\bf{S}}}}}}}}}_{i})}^{2}\): this imposes a constraint on this singlet manifold, energetically penalizing states with a finite value of the total spin S, including states with finite magnetization. As a result, the ground state of the total Hamiltonian is given by the ground state of the short-range Hamiltonian projected on the singlet manifold. This is similar to the analysis in ref. 24, where resonating valence bond (RVB) states with singlets of arbitrarily large size were used as variational wavefunctions. We will denote this state as a long-range singlet state (LRS). Finally, for J = 0, different scenarios are possible: for α large enough, only nearest-neighbouring sites experience an appreciable interaction, and therefore Néel-like order is expected. For smaller values of α, the frustrating nature of the interaction is expected to penalize AFM order, thus favoring disordered phases. This was shown to be the case for α = 3 on the triangular lattice8, and on the square lattice9 (although only for spatially anisotropic interactions in the latter case), where a QSL phase was found.

Summarizing, by varying γ/J and α, we expect three kinds of phases: (i) Néel-like AFM, (ii) a disordered QSL phase, and (iii) an LRS phase. This is substantiated below using a Schwinger–Boson approach, which is capable of capturing all the phases mentioned above. In particular, it provides a natural interpolation scheme between the well-understood RVB and Néel physics discussed above.

In order to unveil the nature of the ground state of the Hamiltonian (1), we apply the bosonic spinon decomposition pioneered in refs. 25,26,27, where the spin operators are represented in terms of new bosonic degrees of freedom, ultimately interpreted as emergent fractional excitations. While this method represents an approximation25,26, it still provides useful information to identify candidate spin liquids. Moreover, the main advantage of this method is its flexibility to interpolate between the different states previously identified. On the one hand, SU(2)-symmetric bosonic ground states are identified with candidate spin liquids. On the other hand, the onset of magnetic order is signalled by the Bose–Einstein condensation of these bosons.

The spin operators on the lattice site j are decomposed as (using sum convention for the Greek indices)

$${{{{{{{{\bf{S}}}}}}}}}_{j}=\frac{1}{2}{b}_{j,\mu }^{{{{\dagger}}} }{{{{{{{{\boldsymbol{\sigma }}}}}}}}}_{\mu \nu }{b}_{j,\nu },$$

where bj,μ is a boson (spinon) with spin μ {↑, ↓}, and σ the vector of Pauli matrices. The mapping is then completed by the constraint \({b}_{j\mu }^{{{{\dagger}}} }{b}_{j\mu }=1\). Insights on the nature of the state are then obtained from the expectation values of the SU(2)-invariant bilinears \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}\,=\,i{\sigma }_{\mu \nu }^{y}\langle {b}_{i\mu }{b}_{j\nu }\rangle /2,\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}=\langle {b}_{i\mu }^{{{{\dagger}}} }{b}_{j\mu }\rangle /2\), which indicate the tendency of the spins at the sites i and j of forming a singlet or to align, respectively. For SU(2)-symmetric states, finite values of \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}\) determine a finite spinon hopping between the lattice sites i and j, thus signalling the emergence of propagating fractional excitations.

After performing a mean-field decoupling of the spinonic Hamiltonian (see ‘Methods’ for further details), the values of \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}\) are self-consistently determined by minimizing the ground-state energy. This task is enabled in practice by using an Ansatz for the values of \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}\). The most natural choice is the manifestly translational-invariant Ansatz \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}={{{{{{{{\mathcal{A}}}}}}}}}_{i-j}\), \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}={{{{{{{{\mathcal{B}}}}}}}}}_{i-j}\), which follows from a projective-symmetry-group analysis28. The resulting saddle-point equations, reported in Eq. (9a, b, c), are reduced to a system of 2N + 1 coupled non-linear equations, for finite-size systems with N = L × L lattice sites. The numerical complexity of the problem still limits the size N of the systems for which a solution can be found.

For finite-size systems, a spontaneous symmetry breaking cannot occur, and therefore the AFM order parameter always vanishes. Accordingly, other criteria are needed to assess the onset of an ordered phase. Here, we identify the onset of an AF-ordered phase when the two following conditions are met: (i) the gap \({E}_{{{{{\mathrm{g}}}}}}\equiv \mathop{\min}\limits_{{{{{{{{\bf{q}}}}}}}}}{E}_{{{{{{{{\bf{q}}}}}}}}}\) in the spinon dispersion closes upon increasing the system size N and (ii) the squared magnetization M2 ≡ ∑jSjS0/N approaches a constant value upon increasing N. Notice that these two indicators also naturally lend themselves to characterize the other phases outlined before: a phase with M2 = 0 corresponds to either a gapped (Eg ≠ 0) or a gapless Eg = 0 QSL, while a phase with M2 ≠ 0 and Eg ≠ 0 can be naturally identified with an LRS state. These criteria are summarized in Table 1.

Table 1 Ground-state phases.

Let us finally discuss the phase diagram in Fig. 1c. The first, main result, is the emergence of a gapped QSL phase (denoted as SL-I) for α 1.25, and γ 5J, characterized by the presence of a gap and by the absence of long-range correlations. This phase appears for any α > 0.05, corresponding to the minimum value here considered, suggesting that the LRS phase is unstable in this region and only exists for α = 0. In addition, our data also show the existence of a gapless QSL phase (denoted SL-II) for intermediate values of α, clearly manifested in the largest available system sizes, as shown in Fig. 2a, b.

Fig. 2: Numerical characterization of the ground state.

a, b Dependence of the spinon gap (upper panel) and square magnetization (lower panel) on the inverse linear system size. The curves refer to values of α and γ/J denoted in Fig. 1c by star symbols, according to the corresponding background colours. The maximum linear size considered is L = 110. Insets: values of the functions \({\xi }_{{E}_{{{{{\mathrm{g}}}}}},{M}^{2}}(L)\) as functions of the linear system size L. c Spinon dispersion for γ = 7J and α = 0.3 (SL-I phase), for given cuts in the first Brillouin zone. Inset: spinon dispersion in the first Brillouin zone. The white lines denote the cuts of the main plot. d Extrapolated gap (blue curve) and square magnetization (red curve) as functions of the exponent α. The background colours reflect the phases illustrated in Fig. 1c. e Spin–spin correlation functions along the lattice axis for different values of the exponent α and for γ = 7J. Inset: spin–spin correlations at short distances.

For γ 5J, the LRS phase is remarkably stable for α 1.25. Here, the system is simultaneously gapped and characterized by long-range correlations (cf. Table 1), which, however, do not correspond to a spontaneous symmetry breaking. Finally, we observe that, as expected, for large values of α, as well as for γ = 0, the system is always in the ordinary Néel–AFM phase.

An example of extrapolated values of M2 and Eg used to build the phase diagram in Fig. 1c is shown in Fig. 2d, as a function of α for γ = 7J. The fitting function used to extrapolate the L →  limit of these observables has the form \({{{{{{{{\mathcal{O}}}}}}}}}_{L}={{{{{{{{\mathcal{O}}}}}}}}}_{\infty }+{b}_{{{{{{{{\mathcal{O}}}}}}}}}{L}^{-{\omega }_{{{{{{{{\mathcal{O}}}}}}}}}}\), with \({{{{{{{{\mathcal{O}}}}}}}}}_{\infty }\), \({b}_{{{{{{{{\mathcal{O}}}}}}}}}\) and \({\omega }_{{{{{{{{\mathcal{O}}}}}}}}}\) fitting parameters. The slightly negative extrapolated values of M2 and Eg are due to the simplified form of the extrapolation function above, which neglects subleading terms in 1/L (cf. ref. 29). This fitting function was identified by a preliminary evaluation of the quantity \({\xi }_{{{{{{{{\mathcal{O}}}}}}}}}(L)=1/{{{{{{\mathrm{ln}}}}}}}\,\left(\frac{{{{{{{{{\mathcal{O}}}}}}}}}_{L-4}-{{{{{{{{\mathcal{O}}}}}}}}}_{L-2}}{{{{{{{{{\mathcal{O}}}}}}}}}_{L-2}-{{{{{{{{\mathcal{O}}}}}}}}}_{L}}\right)\), which displays a linear behaviour in L for algebraic finite-size scaling, while it saturates for an exponential one30. The algebraic finite-size scaling occurring also for gapped phases is imprinted by the algebraic character of the interactions31. For the same reason, the spin–spin correlation functions in the QSL phases also display an algebraically decaying behaviour, rather than the usual short-range one, with an exponent depending continuously on the interaction’s exponent α (cf. Fig. 2e). Algebraic correlations were similarly found for gapped, disordered phases in spin chains with long-range interactions32,33,34, further substantiating the generality of this mechanism.

Besides gap and long-range order, we provide a further observable to characterize the phases here identified, i.e., the dynamical structure factor \({S}_{{{{{{{{\bf{q}}}}}}}}}(\omega )={\int}_{t}{e}^{i\omega t}\langle {{{{{{{{\bf{S}}}}}}}}}_{-{{{{{{{\bf{q}}}}}}}}}(t)\cdot {{{{{{{{\bf{S}}}}}}}}}_{{{{{{{{\bf{q}}}}}}}}}(0)\rangle\), with Sq the Fourier transform of the spin operators with momentum q. Sq(ω), which can be straightforwardly computed from the spinon decomposition35, leads to markedly different features depending on the phase. For the SL-I and SL-II phases (Fig. 3a, b, respectively), the DSF features a broadening originated in the continuum of fractional excitations. On the contrary, the AFM phase (Fig. 3c) shows a sharper signal close to the gapless quasi-particle dispersion, corresponding to the magnonic dispersion expected in the AFM phase. We emphasize that the presence of a gap in the DSF for the SL-II phase is a finite-size effect, and it is expected to close in the thermodynamic limit. Finally, the LRS phase (Fig. 3d) features a broadening similar to the SL-I phase, suggesting the presence of fractionalized excitations.

Fig. 3: Dynamical structure factor.

Sq(ω) as a function of the frequency ω and of the momentum q. Results are shown for: SL-I [panel a, γ = 9J, α = 0.3], SL-II [panel b, γ = 9J, α = 0.6], AFM [panel c, γ = 5J, α = 1.0] and LRS [panel d, γ = 4.4J, α = 0.35].

A final word of caution concerns the accuracy of the bosonic spinon decomposition used here. As a mean-field theory, it provides a qualitatively correct topology of the phase diagram, while the phase borders cannot be expected to be quantitatively accurate.


The Hamiltonian Eq. (1) (or variations of it) can be realized in quantum simulators using trapped ions or ultracold atoms8,36,37. While these platforms provide unprecedented controllability, the realization of low-temperature strongly correlated phases remains challenging. On the converse, solid-state platforms naturally feature strongly correlated physics at cryogenically accessible temperatures. Moreover, the controllability in 2D materials is progressing fast, making them, among others, candidates for quantum simulators38. In the following, we will focus on a scheme for implementing the Hamiltonian (1) in a solid-state system.

Our proposal uses two electronic orbital degrees of freedom, constituting a pseudospin of length S = 1/2. In the absence of a cavity, the pseudospins are assumed to be described by a short-range AFM Heisenberg model, emerging as a strong Mott limit of a Hubbard model for the electronic degrees of freedom: this is the case, e.g., of iridate and ruthenate materials39,40,41. We assume SU(2) symmetry for the sake of simplicity.

As substantiated further below, the coupling of the localized electronic states to the cavity will result in a coupling between the pseudospins and quantized effective magnetic fields. The setup we consider is sketched in Fig. 1a. Two aspects of the long-range Hamiltonian (1) are essential to unveil QSL phases: (i) an AFM character of the induced interaction and (ii) a high degree of symmetry, ideally SU(2). In order to control the symmetry of the emerging cavity-mediated interaction, we propose to use two cavity modes. While a single mode is sufficient to mediate a U(1)-symmetric interaction, a second mode allows for an enhancement to SU(2) symmetry. The required selectivity in the cavity–spin coupling can be achieved via an auxiliary third band which is driven far off-resonance by a laser [see Fig. 1b]. The resulting two-photon transitions involve virtual excitations to the third band and back to one of the two bands implementing the pseudospin degree of freedom. The sign of the cavity-mediated interaction is then finally determined by the detuning between the laser and each cavity mode.

The paramagnetic and diamagnetic coupling terms between electrons and electromagnetic field are given by:

$${H}_{{{{{{{{\rm{int}}}}}}}}}=\frac{1}{2m}{\int}_{{{{{{{{\bf{r}}}}}}}}}{\psi }^{{{{\dagger}}} }({{{{{{{\bf{r}}}}}}}})\left[\right.2e\ {{{{{{{\bf{p}}}}}}}}\cdot {{{{{{{\bf{A}}}}}}}}({{{{{{{\bf{r}}}}}}}},t)+{e}^{2}{{{{{{{{\bf{A}}}}}}}}}^{2}({{{{{{{\bf{r}}}}}}}},t)\left]\right.\psi ({{{{{{{\bf{r}}}}}}}})\ ,$$

with ψ the electronic operators, e the electronic charge and A(r, t) the vector potential. A(r, t) includes the external laser with frequency ωL and the cavity modes a, with frequencies ω,. By choosing the proper polarization for the cavity modes and the laser, the scheme depicted in Fig. 1b can be realized: the laser and the cavity mode a induces transitions between orbitals 1 and 3, while the cavity mode a couples only orbitals 2 and 3. As we assume the electrons to be localized by the strong interaction between particles, due to the strong localization, the field operators can be conveniently expanded onto localized orbitals23: ψ(r) = ∑i,b = 1,2,3wib(r)cib, where wib(r) = wb(r − ri), with ri the position of the centre of the unit cell. Here, the index i runs over the lattice sites and b is the band index. The interaction Hamiltonian (3) thus reads (see ‘Methods’ for further details)

$${H}_{{{\mbox{int}}}}\,=\,\mathop{\sum}\limits_{i}\,\left[{c}_{i3}^{{{{\dagger}}} }{c}_{i1}\left(\right.{\rho }_{31}^{L}+{a}_{\perp }{\rho }_{31}^{\perp }\left)\right.+{c}_{i3}^{{{{\dagger}}} }{c}_{i2}{a}_{\parallel }{\rho }_{32}^{\parallel }+\,{{\mbox{h.c.}}}\,\right],$$

where we neglected counter-rotating terms and changed to the frame rotating with the laser frequency, where \({c}_{i3}\to {c}_{i3}{e}^{-i{\omega }_{L}t}\), \({a}_{{{{{\mathcal{\ell}}}}} }\to {a}_{{{{{\mathcal{\ell}}}}} }{e}^{-i{\omega }_{L}t}\). Correspondingly, the third electron band and the fundamental frequencies of the cavity modes \(\omega_{{{{{\mathcal{\ell}}}}}}\) are shifted as Δ3 = ϵ3 − ωL and \(\Delta_{{{{{\mathcal{\ell}}}}}}\) = \(\omega_{{{{{\mathcal{\ell}}}}}}\) − ωL. The matrix elements \({\rho }_{bb^{\prime} }^{{{{{\mathcal{\ell}}}}} }\), \({{{{{\mathcal{\ell}}}}}}\) {L, , } correspond to the transition rates between the bands b and \(b^{\prime}\).

The effective cavity–spin coupling is then obtained by eliminating the third band adiabatically, assuming the band detuning Δ3 to be much larger than the matrix elements \({\rho }_{bb^{\prime} }^{{{{{\mathcal{\ell}}}}} }\), \({{{{{\mathcal{\ell}}}}}}\) {L, , } and the cavity detunings Δ,. The resulting interaction Hamiltonian describes spins coupled to global, quantum mechanically fluctuating effective magnetic fields:


with \({B}^{z}=-({\rho }_{13}^{L}{\rho }_{31}^{\perp }{a}_{\perp }+\,{{\mbox{h.c.}}}\,)/{{{\Delta }}}_{3}\), \({B}^{x}=-({\rho }_{13}^{L}{\rho }_{32}^{\parallel }{a}_{\parallel }+\,{{\mbox{h.c.}}}\,)/{{{\Delta }}}_{3}\) and \({B}^{y}=-(i{\rho }_{13}^{L}{\rho }_{32}^{\parallel }{a}_{\parallel }+\,{{\mbox{h.c.}}}\,)/{{{\Delta }}}_{3}\), and \({{{{{{{{\bf{S}}}}}}}}}_{i}\,=\,{c}_{ib}^{{{{\dagger}}} }{{{{{{{{\boldsymbol{\sigma }}}}}}}}}_{bb^{\prime} }{c}_{ib^{\prime} }/2\) is the pseudo-spin operator. The values of the effective fluctuating effective magnetic fields Ba, a = x, y, z, reflect the laser-assisted processes illustrated in Fig. 1. For instance, Bx and By, which couple the first and second orbital, result from the laser-assisted excitation of an electron from the first to the third auxiliary band, followed by a decay to the second band with the emission of a cavity photon. The U(1) symmetry of the Hamiltonian results from neglecting the counter-rotating terms, and it is evident from the fact that an excitation from the first to the second band is accompanied only by the creation of a cavity photon, and vice versa. Equation (5) is one of the main results of this paper: the effective quantum magnetic fields Ba couple to all the spins, generating an effective long-range coupling. To further consolidate this insight, we integrate out the cavity field at the level of the Heisenberg equations and obtain an effective Hamiltonian for the spins only42:

$${H}_{{{{{{{{\rm{int}}}}}}}}}=\mathop{\sum}\limits_{ij}\left[{\gamma }_{z}{S}_{i}^{z}{S}_{j}^{z}+{\gamma }_{\perp }({S}_{i}^{x}{S}_{j}^{x}+{S}_{i}^{y}{S}_{j}^{y})\right],$$

with the long-range exchange \({\gamma }_{z}=| {\rho }_{13}^{L}{\rho }_{13}^{\perp }{| }^{2}/({{{\Delta }}}_{3}^{2}{{{\Delta }}}_{\perp })\) and \({\gamma }_{\perp }=| {\rho }_{13}^{L}{\rho }_{23}^{\parallel }{| }^{2}/({{{\Delta }}}_{3}^{2}{{{\Delta }}}_{\parallel })\). The interaction is thus naturally U(1)-symmetric, and full SU(2) symmetry can be achieved by adjusting the cavity-mode detunings. Importantly, by choosing the latter to be positive (i.e. a blue-detuned laser), the cavity-mediated interaction is AFM.

We now briefly show how multi-mode cavities can generate spatially dependent effective spin–spin interactions. To this end, we consider a cavity with a large number of modes. For simplicity, we assume them to correspond to photons propagating as plane waves along the transverse direction with a dispersion \({{{\Delta }}}_{{{{{{\mathcal{\ell}}}}}} ,{{{{{{{\bf{q}}}}}}}}}=\sqrt{{\omega }_{{{{{\mathcal{\ell}}}}} }^{2}+{(c{{{{{{{\bf{q}}}}}}}})}^{2}}-{\omega }_{L}\), with c the speed of light in the medium. The form of the Hamiltonian (5) is then preserved, with the fluctuating magnetic fields now possessing a spatial structure according to

$${B}^{a}=\mathop{\sum}\limits_{{{{{{{{\bf{q}}}}}}}}}{g}_{{{{{{{{\bf{q}}}}}}}}}^{a}\ {a}_{{{{{{{{\bf{q}}}}}}}}}\ {{{{{{{{\rm{e}}}}}}}}}^{i{{{{{{{\bf{q}}}}}}}}\cdot {{{{{{{{\bf{r}}}}}}}}}_{j}}+\,{{\mbox{h.c.}}}\,,$$

with \({g}_{{{{{{{{\bf{q}}}}}}}}}^{a}\) the momentum-dependent version of the coupling reported below Eq. (5). By integrating out the cavity photons, one obtains an effective Hamiltonian as in Eq. (6), where the effective exchange interaction between the spins Sa and Sb is given by \({{{\Gamma }}}_{ij}^{ab}\,=\,{\sum }_{{{{{{{{\bf{q}}}}}}}}}{g}_{{{{{{{{\bf{q}}}}}}}}}^{a}{g}_{{{{{{{{\bf{q}}}}}}}}}^{b}\ {{{{{{{{\rm{e}}}}}}}}}^{-i{{{{{{{\bf{q}}}}}}}}\cdot {{{{{{{{\bf{r}}}}}}}}}_{ij}}/{{{\Delta }}}_{\ell ,{{{{{{{\bf{q}}}}}}}}}.\) While the precise form of \({{{\Gamma }}}_{ij}^{ab}\) depends on the details of \({g}_{{{{{{{{\bf{q}}}}}}}}}^{a}\), its spatial structure is expected to be long-ranged. In fact, the length scale governing the spatial behaviour is proportional to \({\Delta}_{{{{{\mathcal{\ell}}}}}}^{-1/2}\): in THz cavities, the ratio between the lattice size and this length scale is of order 10−4, see, e.g., ref. 14, and, therefore \({{{\Gamma }}}_{ij}^{ab}\), can be effectively modelled as a slowly decaying function. For photonic crystal cavities, the form of \({{{\Gamma }}}_{ij}^{ab}\) can be even further engineered by exploiting the band dispersion of the cavity photons43. The precise form of this function is not expected to qualitatively affect the phase diagram. Accordingly, we choose to parametrize the interaction as \({{{\Gamma }}}_{ij}^{ab}\simeq | {{{{{{{{\bf{r}}}}}}}}}_{i}-{{{{{{{{\bf{r}}}}}}}}}_{j}{| }^{-\alpha }\), with the value of α compactly encoding the interaction range. The values of α achievable with realistic cavity parameters are of order 10−1, and therefore favourable to observe the SL phases (see ‘Methods’ for further details).

We finally provide an estimate for the values of γ in Eq. (1) achievable with this setup (see ‘Methods’ for further details). The dipole matrix elements can be estimated assuming a lattice spacing of few angstroms. For THz cavities with a compression factor of ~10−5 or smaller, a drive with an intensity of ~10 MW cm−2 leads to values of γ of order ~100 K. This number is comparable or larger than typical couplings in antiferromagnets, which range from ~5 K for vanadates44 to ~600 K for iridates45. For α-RuCl3, the (ferromagnetic) Heisenberg interaction is ~40 K, while the Kitaev one is ~80 K, see ref. 46. Accordingly, the spin–liquid phases predicted in the phase diagram in Fig. 1c are achievable with current setups.


In this work, we showed that long-range spin–exchange interactions can be robustly induced by coupling a strongly correlated electron system to the quantum fluctuations of a driven cavity. The electron–cavity coupling gives rise to a variety of tunable spin interactions, including frustrated ones. The thus created cavity-mediated frustration can destroy the magnetic order, favoring disordered spin–liquid states, absent in the cavity-less configuration. We have demonstrated this for an ordinary Heisenberg antiferromagnet, whose ground state manifests an extensive and robust quantum spin–liquid phase when coupled to a cavity. Our results open avenues for engineering quantum spin liquids, sparking the challenge to devise new schemes to control electronic degrees of freedom with quantum light, and to uncover phases of matter that are usually inaccessible. This also represents an exciting perspective for the experimental detection of strongly correlated phases: photons emitted from the cavities carry signatures of the quantum many-body state, which become accessible to standard optical measurements. Our findings are immediately relevant also for quantum simulations. Artificial spin systems with tunable long-range interactions can be currently created using either trapped ions36,47 or ultracold atoms coupled to an optical cavity48,49,50,51. These platforms represent, therefore, ideal candidates to simulate quantum spin liquid phases.


Saddle-point equations for bosonic spinons

In this section, we outline the derivation of the saddle point equations for the spinon bilinear expectation values \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}\).

The spin-exchange terms appearing in Eq. (1) can be recast as \({{{{{{{{\bf{S}}}}}}}}}_{i}{{{{{{{{\bf{S}}}}}}}}}_{j}=:\,{B}_{ij}^{{{{\dagger}}} }{B}_{ij}\,:-{A}_{ij}^{{{{\dagger}}} }{A}_{ij}\) for i ≠ j, where \({A}_{ij}=i{\sigma }_{\mu \nu }^{y}{b}_{i\mu }{b}_{j\nu }/2\) and \({B}_{ij}={b}_{i\mu }^{{{{\dagger}}} }{b}_{j\mu }/2\) are SU(2)-invariant spinonic bilinears. A finite expectation value of these operators indicates the tendency of the spins at the sites i and j of forming a singlet (Aij) or to align (Bij); moreover, it induces a finite bosonic hopping rate between the lattice sites i and j, signalling the existence of propagating fractional excitations. In order to solve for the value of these quantities, we build on the approach of ref. 25. First, the bosonized version of the Hamiltonian (1) is represented as a path integral, with the constraint implemented by a space- and time-dependent Lagrange multiplier λi(t). After decoupling the bilinear products by using a Hubbard–Stratonovich transformation, the expectation values \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}=\langle {A}_{ij}\rangle\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}=\langle {B}_{ij}\rangle\) are obtained as saddle point values of the corresponding action. This approximation imposes the constraint only on average, and the now position- and time-independent Lagrange multiplier λ has to be determined self-consistently. This approximation is equivalent to decoupling the Hamiltonian (1) in bosonic bilinears as:

$$H=\frac{1}{2}\mathop{\sum}\limits_{i,j}\left({\epsilon }_{ij}{b}_{i\mu }^{{{{\dagger}}} }{b}_{j\mu }+i{{{\Delta }}}_{ij}^{* }{\sigma }_{\mu \nu }^{y}{b}_{i\mu }{b}_{j\nu }\right)+\,{{\mbox{h.c.}}}\,+{\varepsilon }_{0},$$

where \({\epsilon }_{ij}={J}_{ij}{{{{{{{{\mathcal{B}}}}}}}}}_{ij}^{* }+{\delta }_{ij}\lambda /2\), \({{{\Delta }}}_{ij}^{* }=-{J}_{ij}{{{{{{{{\mathcal{A}}}}}}}}}_{ij}^{* }\) and \({\varepsilon }_{0}={\sum }_{i,j}(-| {{{{{{{{\mathcal{B}}}}}}}}}_{ij}{| }^{2}+| {{{{{{{{\mathcal{A}}}}}}}}}_{ij}{| }^{2})-2SN\lambda\). As discussed in the main text, we assume a translational-invariant ansatz, i.e. \({{{{{{{{\mathcal{A}}}}}}}}}_{ij}={{{{{{{{\mathcal{A}}}}}}}}}_{i-j}\) and \({{{{{{{{\mathcal{B}}}}}}}}}_{ij}={{{{{{{{\mathcal{B}}}}}}}}}_{i-j}\), able to interpolate between all the expected phases. The two degenerate eigenvalues of H are given by \({E}_{{{{{{{{\bf{q}}}}}}}}}^{2}={\epsilon }_{{{{{{{{\bf{q}}}}}}}}}^{2}-| {{{\Delta }}}_{{{{{{{{\bf{q}}}}}}}}}{| }^{2}\), with ϵq and Δq the Fourier transform of the functions appearing in Eq. (8). By minimizing the ground-state energy \({E}_{0}={\sum }_{{{{{{{{\bf{q}}}}}}}}}\,\left({E}_{{{{{{{{\bf{q}}}}}}}}}-{\epsilon }_{{{{{{{{\bf{q}}}}}}}}}\right)+{\varepsilon }_{0}\) with respect to the variational parameters ϵq, Δq and λ, one obtains the saddle-point equations:

$$1=\frac{1}{2N}\mathop{\sum}\limits_{{{{{{{{\bf{q}}}}}}}}}\frac{{\epsilon }_{{{{{{{{\bf{q}}}}}}}}}}{{E}_{{{{{{{{\bf{q}}}}}}}}}},$$
$${\epsilon }_{{{{{{{{\bf{p}}}}}}}}}=\lambda +\frac{1}{2N}\mathop{\sum}\limits_{{{{{{{{\bf{q}}}}}}}}}{J}_{{{{{{{{\bf{p}}}}}}}}-{{{{{{{\bf{q}}}}}}}}}\left(\frac{{\epsilon }_{{{{{{{{\bf{q}}}}}}}}}}{{E}_{{{{{{{{\bf{q}}}}}}}}}}-1\right),$$
$${{{\Delta }}}_{{{{{{{{\bf{p}}}}}}}}}=\frac{1}{2N}\mathop{\sum}\limits_{{{{{{{{\bf{q}}}}}}}}}{J}_{{{{{{{{\bf{p}}}}}}}}-{{{{{{{\bf{q}}}}}}}}}\frac{{{{\Delta }}}_{{{{{{{{\bf{q}}}}}}}}}}{{E}_{{{{{{{{\bf{q}}}}}}}}}},$$

with Jq the Fourier transform of Jij. These equations provide the full momentum dependence of the functions ϵq and Δq. The actual number of unknowns increases with the range of the interaction. In fact, for short-range interactions, the momentum dependence can be found analytically, and only a few parameters are left to be computed self-consistently. For long-range interactions, instead, the full momentum dependence needs to be found numerically. Equation (9a, b, c) amounts to a system of 2N + 1 coupled non-linear equations, with N the total number of sites. To find the roots of these equations, we used a trust region solver as provided by the Julia NLSolve library, with the accepted residual norm set to 10−8. The error of the numerical solution for a finite system of N sites is, therefore, negligible compared to the extrapolation to the thermodynamic limit. Determining \({\xi }_{{{{{{{{\mathcal{O}}}}}}}}}(L)\) by a least-squares fit resulted in relative errors between 0.5 and 2%, which we consider to be sufficient for the analysis performed here. The root-finding algorithm is accelerated by exploiting vectorization for the evaluation of the saddle-point equations where possible, and by parallelization via OpenBLAS.

Implementation details

Here, we provide additional details to the set-up described in the main text. The vector potential can be written as

$${{{{{{{\bf{A}}}}}}}}({{{{{{{\bf{r}}}}}}}},t)={{{\Omega }}}_{L}{{{{{{{{\bf{u}}}}}}}}}_{L}{\varphi }_{L}({{{{{{{\bf{r}}}}}}}}){e}^{i{\omega }_{L}t}+\mathop{\sum}\limits_{{{{{{\mathcal{\ell}}}}}} =\parallel ,\perp }{{{{{{{{\mathcal{N}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }{{{{{{{{\bf{u}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }{a}_{{{{{\mathcal{\ell}}}}} }{\varphi }_{{{{{\mathcal{\ell}}}}} }({{{{{{{\bf{r}}}}}}}})+{{{{{{{\rm{h.c.}}}}}}}},$$

where \({{{\Omega }}}_{L}^{2}\) and ωL denote the laser intensity and frequency, respectively. Here u and φ(r) are the polarization vector and the mode wavefunction. For the cavity modes, labelled by \({{{{{\mathcal{\ell}}}}}}\) = , , the wave function is normalized over the finite volume Vc and \({{{{{{{{\mathcal{N}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }=\sqrt{1/2{\omega }_{{{{{\mathcal{\ell}}}}} }{\epsilon }_{0}{\epsilon }_{r}}\), where \(\omega_{{{{{\mathcal{\ell}}}}}}\) is the mode fundamental frequency and ϵ0, ϵr are the vacuum and relative permittivity of the material, respectively.

We assume that that the mode wavefunctions φ(r) does not vary significantly over the extent of the Wannier functions. By tuning the polarization vectors u to selectively couple the orbitals as in Fig. 1, and by performing the rotating-wave approximation, the resulting paramagnetic Hamiltonian term is given by Eq. (4), with

$${\rho }_{b{b}^{\prime}}^{{{{{\mathcal{\ell}}}}} }=\frac{e{{{{{{{{\mathcal{N}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }}{m}{\varphi }_{{{{{\mathcal{\ell}}}}} }{{{{{{{{\bf{u}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }\cdot \langle {w}_{ib}| {{{{{{{\bf{p}}}}}}}}| {w}_{ib^{\prime} }\rangle .$$

The expression for \({\rho }_{b{b}^{\prime}}^{L}\) can be obtained from the previous equation by replacing \({{{{{{{{\mathcal{N}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }\) with ΩL.

The diamagnetic part of the Hamiltonian (3) reads, after neglecting higher-order electron-photon processes of the type \({c}_{3}^{{{{\dagger}}} }{c}_{3}(a+{a}^{{{{\dagger}}} })\) and \({c}_{3}^{{{{\dagger}}} }{c}_{3}{a}^{{{{\dagger}}} }a\):

$${H}_{{{\mbox{int,dia}}}}=\mathop{\sum}\limits_{{{{{{\mathcal{\ell}}}}}} =\perp ,\parallel }{\delta }_{{{{{\mathcal{\ell}}}}} }\ {a}_{{{{{\mathcal{\ell}}}}} }^{{{{\dagger}}} }{a}_{{{{{\mathcal{\ell}}}}} }+\mathop{\sum}\limits_{i}{\delta }_{3}\ {c}_{i3}^{{{{\dagger}}} }{c}_{i3},$$

plus a term linear in the cavity fields, which vanishes as the laser and cavity wavefunctions are orthogonal. The shifts

$${\delta }_{{{{{\mathcal{\ell}}}}} }=\frac{{e}^{2}}{m}{V}_{{{\mbox{e}}}}\ {{{{{{{{\mathcal{N}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }^{2}| {\varphi }_{{{{{\mathcal{\ell}}}}} }{| }^{2},$$
$${\delta }_{3}=\frac{{e}^{2}{{{\Omega }}}_{L}^{2}}{m}| {\varphi }_{L}{| }^{2}+\mathop{\sum}\limits_{{{{{{\mathcal{\ell}}}}}} =\perp ,\parallel }\frac{{e}^{2}{{{{{{{{\mathcal{N}}}}}}}}}_{{{{{\mathcal{\ell}}}}} }^{2}}{2m}| {\varphi }_{{{{{\mathcal{\ell}}}}} }{| }^{2},$$

renormalize the energies of the cavity modes and of the third band, respectively. By assuming that the band detuning \(| {\widetilde{{{\Delta }}}}_{3}|\) is much larger than the coupling strengths and the cavity detunings \({\widetilde{{{\Delta }}}}_{{{{{\mathcal{\ell}}}}} }\), the third band can be adiabatically eliminated, leading to Eq. (5) in the main text, including an additional term \({B}_{0}\;{\sum }_{i}{S}_{i}^{z}\), with \({B}_{0}=| {\rho }_{13}^{L}{| }^{2}/{\widetilde{{{\Delta }}}}_{3}\). This effective classical magnetic field breaks explicitly the SU(2) symmetry, but it is much smaller than the spin exchange and therefore it can be safely neglected.

Estimate of interaction strength and range

We consider a THz laser (ωL = 100 THz) with intensity \({{{\Omega }}}_{L}^{2}=10\) MW cm−2, with a small detuning from the cavity frequency Δ = Δ = 10−2 THz. The compression factor of the cavity is assumed to be Λ = 10−5. The detuning from the third band is Δ3 = 1 THz, thus satisfying the condition Δ3 Δ. We estimate the matrix as follows: \(\langle {w}_{i1}| {{{{{{{\bf{p}}}}}}}}| {w}_{i3^{\prime} }\rangle \sim m{\omega }_{13}\langle {w}_{i1}| {{{{{{{\bf{r}}}}}}}}| {w}_{i3^{\prime} }\rangle\), with ω13 = ωL + ΔL and \(\langle {w}_{i1}| {{{{{{{\bf{r}}}}}}}}| {w}_{i3^{\prime} }\rangle =10\,\mathring{\rm A}\), the same order of magnitude of a typical lattice spacing. Using the formulas derived in the text, one then estimates a long-range interaction with strength γ ~ 100 K.

We also provide an estimate of the values of α. To this end, we evaluate the explicit form of Γ(rij) as reported in the text below Eq. (7). For simplicity, we assume Δ = Δ ≡ Δ and ω = ω ≡ ω. The Rabi-like couplings gq inherits the momentum dependence from the normalization of every mode, i.e. \({g}_{{{{{{{{\bf{q}}}}}}}}}\propto {({\omega }^{2}+{(c{{{{{{{\bf{q}}}}}}}})}^{2})}^{-1/4}\). Accordingly, the cavity-mediated exchange is given by:

$${{\Gamma }}({{{{{{{{\bf{r}}}}}}}}}_{ij})\propto \mathop{\sum}\limits_{{{{{{{{\bf{q}}}}}}}}}\frac{{{{{{{{{\rm{e}}}}}}}}}^{-i{{{{{{{\bf{q}}}}}}}}\cdot {{{{{{{{\bf{r}}}}}}}}}_{ij}}}{\sqrt{{\omega }^{2}+{(c{{{{{{{\bf{q}}}}}}}})}^{2}}(\sqrt{{\omega }_{c}^{2}+{(c{{{{{{{\bf{q}}}}}}}})}^{2}}-{\omega }_{L})}.$$

The corresponding integral is computed numerically, and the results shown in Fig. 4 over a range of 50 lattice sites, for different values of the cavity detuning. The values of α obtained are reported in the figure.

Fig. 4: Normalized long-range spin-exchange Γ(r) as a function of the distance between sites r.

The cavity fundamental wavelength is λ = 2πc/ω = 105a, with a the lattice spacing. The symbols correspond to different values of the cavity detuning: Δ = 0.5ω (orange), Δ = 0.1ω (red), Δ = 0.05ω (purple) and Δ = 0.01ω (green). The dashed curves correspond to the fitted power laws, with exponents reported in the legend.

Heating effects

A possible advantage of our scheme is that it does not rely on the laser being resonant with any electronic or phononic excitation, with the only tradeoff of a decreasing coupling strength as the detuning increases. Exploiting the variability of the detuning, as well as the knowledge of the relevant excitation modes of the quantum material and the cavity, resonances can be avoided and heating is pushed to later times. This said, we can on the other hand consider the worst-case scenario, where heating does efficiently take place, and estimate the amount of it, following ref. 52. As a paradigmatic material, we consider RuCl3, which features orbital pseudospin, and whose low-temperature properties have been intensely studied46. In order to estimate the heating of the material surface due to the external laser, we evaluate the energy density deposited by the laser using the formula:

$$\epsilon (z)=(1-R)\frac{{{{{{{{\mathcal{F}}}}}}}}}{{d}_{{{\mbox{p}}}}}{{{{{{{{\rm{e}}}}}}}}}^{-\frac{z}{{d}_{{{\mbox{p}}}}}},$$

with z the depth in the sample, R ~ 0.05 the material reflectivity at the laser frequency considered, i.e. 100 THz (cf. ref. 53), dp the penetration depth, chosen of an order of micrometres and \({{{{{{{\mathcal{F}}}}}}}}\) is the excitation energy density. In order to achieve lasing in the desired frequency range, a short-pulse protocol can be used. By considering pulses of ~10 ps, and maximal laser intensity of 10 MW cm−2, we require \({{{{{{{\mathcal{F}}}}}}}}=1{0}^{-5}\) J cm−2. In order to estimate the increase in temperature due to the deposited energy density ϵ(z), we assume that thermalization time is fast, and use the following relation:

$$\epsilon (z)=\frac{1}{{V}_{{{{{\mathrm{m}}}}}}}\int\limits_{{T}_{0}}^{{T}_{{{{{\mathrm{f}}}}}}(z)}{C}_{{{{{\mathrm{p}}}}}}(T){{{{{{{\rm{d}}}}}}}}T,$$

where Vm = 53.32 cm3 mol−1 is the molar volume of α-RuCl3, and Cp is the molar heat capacity, for which we use the value fitted from the measurements of ref. 54, and T0 is the initial temperature in the sample. For an initial temperature of T0 = 2 K, the rise in temperature as a function of z and for different values of the penetration length dp are reported in Fig. 5 here. The estimated temperature increase in the first layers is in between 5 and 20 K, depending on the penetration depth. In order to understand the impact of heating on the candidate QSL phases, a simple criterion compares the temperature increase with the gap (at least for the gapped spin liquid phase, which we dubbed SL-I in the main text). Robustness of the QSL phase then requires the temperature increase to be smaller than the gap, whose scale is given by the material couplings and by the cavity-induced interaction. For α-RuCl3, as discussed in the main text, interactions lie between 40 and 80 K. Accordingly, the heating induced by the laser is not expected to destabilize the gapped QSL phase. For what concerns the gapless QSL phase (SL-II in the main text), the previous argument is clearly inapplicable, and its robustness against thermal fluctuations must be assessed using a more sophisticated approach, e.g. by solving the saddle-point equations at finite temperature.

Fig. 5: Heating effects.

Laser-induced temperature increase as a function of the depth in the sample, for different values of the laser penetration depths.

Code availability

The code that supports the plots within this paper is available from the corresponding author upon reasonable request.


  1. 1.

    Savary, L. & Balents, L. Quantum spin liquids: a review. Rep. Prog. Phys. 80, 016502 (2016).

    ADS  PubMed  Article  CAS  Google Scholar 

  2. 2.

    Zhou, Y. I., Kanoda, K. & Ng, Tai-Kai Quantum spin liquid states. Rev. Mod. Phys. 89, 025003 (2017).

    ADS  MathSciNet  Article  Google Scholar 

  3. 3.

    Knolle, J. & Moessner, R. A field guide to spin liquids. Annu. Rev. Condens. Matter Phys. 10, 451–472 (2019).

    ADS  Article  CAS  Google Scholar 

  4. 4.

    Balents, L. Spin liquids in frustrated magnets. Nature 464, 199–208 (2010).

    ADS  PubMed  Article  CAS  Google Scholar 

  5. 5.

    Castelnovo, C., Moessner, R. & Sondhi, S. L. Spin ice, fractionalization, and topological order. Annu. Rev. Condens. Matter Phys. 3, 35–55 (2012).

    Article  CAS  Google Scholar 

  6. 6.

    Gingras, M. J. P. & McClarty, P. A. Quantum spin ice: a search for gapless quantum spin liquids in pyrochlore magnets. Rep. Prog. Phys. 77, 056501 (2014).

    ADS  PubMed  Article  CAS  Google Scholar 

  7. 7.

    Norman, M. R. Colloquium: Herbertsmithite and the search for the quantum spin liquid. Rev. Mod. Phys. 88, 041002 (2016).

    ADS  MathSciNet  Article  Google Scholar 

  8. 8.

    Yao, N. Y., Zaletel, M. P., Stamper-Kurn, D. M. & Vishwanath, A. A quantum dipolar spin liquid. Nat. Phys. 14, 405–410 (2018).

    Article  CAS  Google Scholar 

  9. 9.

    Zou, H., Zhao, E. & Liu, W. V. Frustrated magnetism of dipolar molecules on a square optical lattice: prediction of a quantum paramagnetic ground state. Phys. Rev. Lett. 119, 050401 (2017).

    ADS  PubMed  Article  Google Scholar 

  10. 10.

    Keleş, A. & Zhao, E. Absence of long-range order in a triangular spin system with dipolar interactions. Phys. Rev. Lett. 120, 187202 (2018).

    ADS  PubMed  Article  Google Scholar 

  11. 11.

    Sentef, M. A., Ruggenthaler, M. & Rubio, A. Cavity quantum-electrodynamical polaritonically enhanced electron-phonon coupling and its influence on superconductivity. Sci. Adv. 4, eaau6969 (2018).

    ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

  12. 12.

    Mazza, G. & Georges, A. Superradiant quantum materials. Phys. Rev. Lett. 122, 017401 (2019).

    ADS  PubMed  Article  CAS  Google Scholar 

  13. 13.

    Andolina, G. M., Pellegrino, F. M. D., Giovannetti, V., MacDonald, A. H. & Polini, M. Cavity quantum electrodynamics of strongly correlated electron systems: a no-go theorem for photon condensation. Phys. Rev. B 100, 121109 (2019).

    ADS  Article  CAS  Google Scholar 

  14. 14.

    Schlawin, F., Cavalleri, A. & Jaksch, D. Cavity-mediated electron-photon superconductivity. Phys. Rev. Lett. 122, 133602 (2019).

    ADS  PubMed  Article  CAS  Google Scholar 

  15. 15.

    Curtis, J. B., Raines, Z. M., Allocca, A. A., Hafezi, M. & Galitski, V. M. Cavity quantum eliashberg enhancement of superconductivity. Phys. Rev. Lett. 122, 167002 (2019).

    ADS  PubMed  Article  CAS  Google Scholar 

  16. 16.

    Thomas, A. et al. Exploring superconductivity under strong coupling with the vacuum electromagnetic field. Preprint at arXiv (2019).

  17. 17.

    Gao, H., Schlawin, F., Buzzi, M., Cavalleri, A. & Jaksch, D. Photoinduced electron pairing in a driven cavity. Phys. Rev. Lett. 125, 053602 (2020).

    ADS  PubMed  Article  CAS  Google Scholar 

  18. 18.

    Ashida, Y. et al. Quantum electrodynamic control of matter: cavity-enhanced ferroelectric phase transition. Phys. Rev. X 10, 041027 (2020).

    CAS  Google Scholar 

  19. 19.

    Schuler, M., Bernardis, DanieleDe, Läuchli, A. M. & Rabl, P. The vacua of dipolar cavity quantum electrodynamics. SciPost Phys. 9, 66 (2020).

    ADS  MathSciNet  Article  Google Scholar 

  20. 20.

    Chakraborty, A. & Piazza, F. Non-bcs-type enhancement of superconductivity from long-range photon fluctuations. Preprint at arXiv (2020).

  21. 21.

    Sentef, M. A., Li, J., Künzel, F. & Eckstein, M. Quantum to classical crossover of floquet engineering in correlated quantum systems. Phys. Rev. Res. 2, 033033 (2020).

    Article  CAS  Google Scholar 

  22. 22.

    Kiffner, M., Coulthard, J. R., Schlawin, F., Ardavan, A. & Jaksch, D. Manipulating quantum materials with quantum light. Phys. Rev. B 99, 085116 (2019).

    ADS  Article  CAS  Google Scholar 

  23. 23.

    Auerbach, A. Interacting Electrons and Quantum Magnetism (Springer Science & Business Media, 2012)

  24. 24.

    Liang, S., Doucot, B. & Anderson, P. W. Some new variational resonating-valence-bond-type wave functions for the spin- antiferromagnetic heisenberg model on a square lattice. Phys. Rev. Lett. 61, 365–368 (1988).

    ADS  PubMed  Article  CAS  Google Scholar 

  25. 25.

    Arovas, D. P. & Auerbach, A. Functional integral theories of low-dimensional quantum heisenberg models. Phys. Rev. B 38, 316–332 (1988).

    ADS  Article  CAS  Google Scholar 

  26. 26.

    Read, N. & Sachdev, S. Large-n expansion for frustrated quantum antiferromagnets. Phys. Rev. Lett. 66, 1773–1776 (1991).

    ADS  PubMed  Article  CAS  Google Scholar 

  27. 27.

    Ceccatto, H. A., Gazza, C. J. & Trumper, A. E. Nonclassical disordered phase in the strong quantum limit of frustrated antiferromagnets. Phys. Rev. B 47, 12329–12332 (1993).

    ADS  Article  CAS  Google Scholar 

  28. 28.

    Wen, Xiao-Gang Quantum orders and symmetric spin liquids. Phys. Rev. B 65, 165113 (2002).

    ADS  Article  CAS  Google Scholar 

  29. 29.

    Sorella, S., Otsuka, Y. & Yunoki, S. Absence of a spin liquid phase in the hubbard model on the honeycomb lattice. Sci. Rep. 2, 992 (2012).

    ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

  30. 30.

    Golinelli, O., Jolicoeur, T. H. & Lacaze, R. Finite-lattice extrapolations for a haldane-gap antiferromagnet. Phys. Rev. B 50, 3037–3044 (1994).

    ADS  Article  CAS  Google Scholar 

  31. 31.

    Campa, A., Dauxois, T., Fanelli, D., & Ruffo, S. Physics of Long-Range Interacting Systems (Oxford Univ. Press, 2014).

  32. 32.

    Hauke, P. et al. Complete devil’s staircase and crystal-superfluid transitions in a dipolarXXZspin chain: a trapped ion quantum simulation. N. J. Phys. 12, 113037 (2010).

    Article  Google Scholar 

  33. 33.

    Koffel, T., Lewenstein, M. & Tagliacozzo, L. Entanglement entropy for the long-range ising chain in a transverse field. Phys. Rev. Lett. 109, 267203 (2012).

    ADS  PubMed  Article  CAS  Google Scholar 

  34. 34.

    Peter, D., Müller, S., Wessel, S. & Büchler, H. P. Anomalous behavior of spin systems with dipolar interactions. Phys. Rev. Lett. 109, 025303 (2012).

    ADS  PubMed  Article  CAS  Google Scholar 

  35. 35.

    Messio, L., Cépas, O. & Lhuillier, C. Schwinger-boson approach to the kagome antiferromagnet with dzyaloshinskii-moriya interactions: phase diagram and dynamical structure factors. Phys. Rev. B 81, 064428 (2010).

    ADS  Article  CAS  Google Scholar 

  36. 36.

    Richerme, P. et al. Non-local propagation of correlations in quantum systems with long-range interactions. Nature 511, 198–201 (2014).

    ADS  PubMed  Article  CAS  Google Scholar 

  37. 37.

    Mivehvar, F., Piazza, F., Donner, T., & Ritsch, H. Cavity qed with quantum gases: new paradigms in many-body physics. Preprint at arXiv (2021).

  38. 38.

    Kennes, D. M. et al. Moiréheterostructures as a condensed-matter quantum simulator. Nat. Phys. 17, 155–163 (2021).

    Article  CAS  Google Scholar 

  39. 39.

    Jackeli, G. & Khaliullin, G. Mott insulators in the strong spin-orbit coupling limit: From heisenberg to a quantum compass and kitaev models. Phys. Rev. Lett. 102, 017205 (2009).

    ADS  PubMed  Article  CAS  Google Scholar 

  40. 40.

    Georges, A., de’ Medici, L. & Mravlje, J. Strong correlations from hund’s coupling. Annu. Rev. Condens. Matter Phys. 4, 137–178 (2013).

    ADS  Article  CAS  Google Scholar 

  41. 41.

    Rau, J. G., Lee, EricKin-Ho & Kee, Hae-Young Spin-orbit physics giving rise to novel phases in correlated systems: Iridates and related materials. Annu. Rev. Condens. Matter Phys. 7, 195–221 (2016).

    ADS  Article  CAS  Google Scholar 

  42. 42.

    Mivehvar, F., Ritsch, H. & Piazza, F. Cavity-quantum-electrodynamical toolbox for quantum magnetism. Phys. Rev. Lett. 122, 113603 (2019).

    ADS  PubMed  Article  CAS  Google Scholar 

  43. 43.

    Joannopoulos, J. D., Johnson, S. G., Winn, J. N., & Meade, R. D. Photonic Crystals: Molding the Flow of Light 2nd edn (Princeton Univ. Press, 2008).

  44. 44.

    Melzi, R. et al. li2VO(si, ge)O4, a prototype of a two-dimensional frustrated quantum heisenberg antiferromagnet. Phys. Rev. Lett. 85, 1318–1321 (2000).

    ADS  PubMed  Article  CAS  Google Scholar 

  45. 45.

    Kim, J. et al. Magnetic excitation spectra of sr2iro4 probed by resonant inelastic x-ray scattering: establishing links to cuprate superconductors. Phys. Rev. Lett. 108, 177003 (2012).

    ADS  PubMed  Article  CAS  Google Scholar 

  46. 46.

    Banerjee, A. et al. Proximate kitaev quantum spin liquid behaviour in a honeycomb magnet. Nat. Mater. 15, 733–740 (2016).

    ADS  PubMed  Article  CAS  Google Scholar 

  47. 47.

    Jurcevic, P. et al. Quasiparticle engineering and entanglement propagation in a quantum many-body system. Nature 511, 202–205 (2014).

    ADS  PubMed  Article  CAS  Google Scholar 

  48. 48.

    Landini, M., Dogra, N., Kroeger, K., Hruby, L., Donner, T. & Esslinger, T. Formation of a spin texture in a quantum gas coupled to a cavity. Phys. Rev. Lett. 120, 223602 (2018).

    ADS  PubMed  Article  CAS  Google Scholar 

  49. 49.

    Kroeze, R. M., Guo, Y., Vaidya, V. D., Keeling, J. & Lev, B. L. Spinor self-ordering of a quantum gas in a cavity. Phys. Rev. Lett. 121, 163601 (2018).

    ADS  PubMed  Article  CAS  Google Scholar 

  50. 50.

    Krešić, I. et al. Spontaneous light-mediated magnetism in cold atoms. Commun. Phys. 1, 1–9 (2018).

    Article  CAS  Google Scholar 

  51. 51.

    Davis, E. J., Bentsen, G., Homeier, L., Li, T. & Schleier-Smith, M. H. Photon-mediated spin-exchange dynamics of spin-1 atoms. Phys. Rev. Lett. 122, 010405 (2019).

    ADS  PubMed  Article  CAS  Google Scholar 

  52. 52.

    Gao, H., Schlawin, F., Buzzi, M., Cavalleri, A. & Jaksch, D. Photoinduced electron pairing in a driven cavity. Phys. Rev. Lett. 125, 053602 (2020).

    ADS  PubMed  Article  CAS  Google Scholar 

  53. 53.

    Reschke, S., Mayr, F., Wang, Z., Do, Seung-Hwan, Choi, K.-Y. & Loidl, A. Electronic and phonon excitations in α − RuCl3. Phys. Rev. B 96, 165120 (2017).

    ADS  Article  Google Scholar 

  54. 54.

    Cao, H. B. et al. Low-temperature crystal and magnetic structure of α − rucl3. Phys. Rev. B 93, 134423 (2016).

    ADS  Article  CAS  Google Scholar 

Download references


We warmly acknowledge discussions with C. Hickey, G. Mazza, A. Rosch, M. Scherer and especially S. Trebst. We acknowledge support by the funding from the European Research Council (ERC) under the Horizon 2020 research and innovation programme, Grant Agreement No. 647434 (DOQS), and by the DFG Collaborative Research Centre (CRC) 1238 Project No. 277146847—projects C02, C03 and C04.


Open Access funding enabled and organized by Projekt DEAL.

Author information




S.D. and A.C. designed the research, D.K. performed the numerical simulations, A.C. and F.P. developed the implementation scheme. C.P.Z. computed and analysed the dynamical structure factor data. All authors analysed the results and contributed to the manuscript.

Corresponding author

Correspondence to Alessio Chiocchetta.

Ethics declarations

Competing interests

The authors declare no competing interests.

Additional information

Peer review informationNature Communications thanks the anonymous reviewers for their contribution to the peer review of this work.

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

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Chiocchetta, A., Kiese, D., Zelle, C.P. et al. Cavity-induced quantum spin liquids. Nat Commun 12, 5901 (2021).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI:


By submitting a comment you agree to abide by our Terms and Community Guidelines. If you find something abusive or that does not comply with our terms or guidelines please flag it as inappropriate.


Quick links

Nature Briefing

Sign up for the Nature Briefing newsletter — what matters in science, free to your inbox daily.

Get the most important science stories of the day, free in your inbox. Sign up for Nature Briefing