## Abstract

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.

## Introduction

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 distances^{1,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 lattices^{4,5,6,7}. Another one proceeds via the energetic competition of couplings of different ranges, like in the antiferromagnetic (AFM) *J*_{1}–*J*_{2} Heisenberg model or dipolar-interacting systems^{8,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.

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 phases^{11,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 way^{12,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.

## Results

### Model

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

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 **r**_{ij} ≡ **r**_{i} − **r**_{j}. 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 order^{23}. 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 lattice^{8}, and on the square lattice^{9} (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 approximation^{25,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)

where *b*_{j,μ} 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 analysis^{28}. The resulting saddle-point equations, reported in Eq. (9a, b, c), are reduced to a system of 2*N* + 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 *M*^{2} ≡ ∑_{j}∣**S**_{j} ⋅ **S**_{0}∣/*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 *M*^{2} = 0 corresponds to either a gapped (*E*_{g} ≠ 0) or a gapless *E*_{g} = 0 QSL, while a phase with *M*^{2} ≠ 0 and *E*_{g} ≠ 0 can be naturally identified with an LRS state. These criteria are summarized in Table 1.

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 *γ* ≳ 5*J*, 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.

For *γ* ≲ 5*J*, 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 *M*^{2} and *E*_{g} used to build the phase diagram in Fig. 1c is shown in Fig. 2d, as a function of *α* for *γ* = 7*J*. 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 *M*^{2} and *E*_{g} 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 one^{30}. The algebraic finite-size scaling occurring also for gapped phases is imprinted by the algebraic character of the interactions^{31}. 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 interactions^{32,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 **S**_{q} the Fourier transform of the spin operators with momentum **q**. *S*_{q}(*ω*), which can be straightforwardly computed from the spinon decomposition^{35}, 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.

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.

### Implementation

The Hamiltonian Eq. (1) (or variations of it) can be realized in quantum simulators using trapped ions or ultracold atoms^{8,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 simulators^{38}. 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 materials^{39,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:

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 orbitals^{23}: *ψ*(**r**) = ∑_{i,b = 1,2,3}*w*_{ib}(**r**)*c*_{ib}, where *w*_{ib}(**r**) = *w*_{b}(**r** − **r**_{i}), with **r**_{i} 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)

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 *B*^{a}, *a* = *x*, *y*, *z*, reflect the laser-assisted processes illustrated in Fig. 1. For instance, *B*^{x} and *B*^{y}, 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 *B*^{a} 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 only^{42}:

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

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 *S*^{a} and *S*^{b} 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 photons^{43}. 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 vanadates^{44} to ~600 K for iridates^{45}. For *α*-RuCl_{3}, 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.

## Discussion

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 ions^{36,47} or ultracold atoms coupled to an optical cavity^{48,49,50,51}. These platforms represent, therefore, ideal candidates to simulate quantum spin liquid phases.

## Methods

### 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 (*A*_{ij}) or to align (*B*_{ij}); 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:

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:

with *J*_{q} the Fourier transform of *J*_{ij}. 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 2*N* + 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

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 *V*_{c} 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

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\):

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

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 Γ(**r**_{ij}) as reported in the text below Eq. (7). For simplicity, we assume Δ_{⊥} = Δ_{∥} ≡ Δ and *ω*_{⊥} = *ω*_{∥} ≡ *ω*. The Rabi-like couplings *g*_{q} 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:

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.

### 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 RuCl_{3}, which features orbital pseudospin, and whose low-temperature properties have been intensely studied^{46}. 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:

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}), *d*_{p} 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:

where *V*_{m} = 53.32 cm^{3} mol^{−1} is the molar volume of *α*-RuCl_{3}, and *C*_{p} is the molar heat capacity, for which we use the value fitted from the measurements of ref. ^{54}, and *T*_{0} is the initial temperature in the sample. For an initial temperature of *T*_{0} = 2 K, the rise in temperature as a function of *z* and for different values of the penetration length *d*_{p} 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 *α*-RuCl_{3}, 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.

## Code availability

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

## References

- 1.
Savary, L. & Balents, L. Quantum spin liquids: a review.

*Rep. Prog. Phys.***80**, 016502 (2016). - 2.
Zhou, Y. I., Kanoda, K. & Ng, Tai-Kai Quantum spin liquid states.

*Rev. Mod. Phys.***89**, 025003 (2017). - 3.
Knolle, J. & Moessner, R. A field guide to spin liquids.

*Annu. Rev. Condens. Matter Phys.***10**, 451–472 (2019). - 4.
Balents, L. Spin liquids in frustrated magnets.

*Nature***464**, 199–208 (2010). - 5.
Castelnovo, C., Moessner, R. & Sondhi, S. L. Spin ice, fractionalization, and topological order.

*Annu. Rev. Condens. Matter Phys.***3**, 35–55 (2012). - 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). - 7.
Norman, M. R. Colloquium: Herbertsmithite and the search for the quantum spin liquid.

*Rev. Mod. Phys.***88**, 041002 (2016). - 8.
Yao, N. Y., Zaletel, M. P., Stamper-Kurn, D. M. & Vishwanath, A. A quantum dipolar spin liquid.

*Nat. Phys.***14**, 405–410 (2018). - 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). - 10.
Keleş, A. & Zhao, E. Absence of long-range order in a triangular spin system with dipolar interactions.

*Phys. Rev. Lett.***120**, 187202 (2018). - 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). - 12.
Mazza, G. & Georges, A. Superradiant quantum materials.

*Phys. Rev. Lett.***122**, 017401 (2019). - 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). - 14.
Schlawin, F., Cavalleri, A. & Jaksch, D. Cavity-mediated electron-photon superconductivity.

*Phys. Rev. Lett.***122**, 133602 (2019). - 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). - 16.
Thomas, A. et al. Exploring superconductivity under strong coupling with the vacuum electromagnetic field. Preprint at

*arXiv*https://arxiv.org/abs/1911.01459 (2019). - 17.
Gao, H., Schlawin, F., Buzzi, M., Cavalleri, A. & Jaksch, D. Photoinduced electron pairing in a driven cavity.

*Phys. Rev. Lett.***125**, 053602 (2020). - 18.
Ashida, Y. et al. Quantum electrodynamic control of matter: cavity-enhanced ferroelectric phase transition.

*Phys. Rev. X***10**, 041027 (2020). - 19.
Schuler, M., Bernardis, DanieleDe, Läuchli, A. M. & Rabl, P. The vacua of dipolar cavity quantum electrodynamics.

*SciPost Phys.***9**, 66 (2020). - 20.
Chakraborty, A. & Piazza, F. Non-bcs-type enhancement of superconductivity from long-range photon fluctuations. Preprint at

*arXiv*https://arxiv.org/abs/2008.06513 (2020). - 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). - 22.
Kiffner, M., Coulthard, J. R., Schlawin, F., Ardavan, A. & Jaksch, D. Manipulating quantum materials with quantum light.

*Phys. Rev. B***99**, 085116 (2019). - 23.
Auerbach, A.

*Interacting Electrons and Quantum Magnetism*(Springer Science & Business Media, 2012) - 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). - 25.
Arovas, D. P. & Auerbach, A. Functional integral theories of low-dimensional quantum heisenberg models.

*Phys. Rev. B***38**, 316–332 (1988). - 26.
Read, N. & Sachdev, S. Large-n expansion for frustrated quantum antiferromagnets.

*Phys. Rev. Lett.***66**, 1773–1776 (1991). - 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). - 28.
Wen, Xiao-Gang Quantum orders and symmetric spin liquids.

*Phys. Rev. B***65**, 165113 (2002). - 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). - 30.
Golinelli, O., Jolicoeur, T. H. & Lacaze, R. Finite-lattice extrapolations for a haldane-gap antiferromagnet.

*Phys. Rev. B***50**, 3037–3044 (1994). - 31.
Campa, A., Dauxois, T., Fanelli, D., & Ruffo, S.

*Physics of Long-Range Interacting Systems*(Oxford Univ. Press, 2014). - 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). - 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). - 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). - 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). - 36.
Richerme, P. et al. Non-local propagation of correlations in quantum systems with long-range interactions.

*Nature***511**, 198–201 (2014). - 37.
Mivehvar, F., Piazza, F., Donner, T., & Ritsch, H. Cavity qed with quantum gases: new paradigms in many-body physics. Preprint at

*arXiv*https://doi.org/10.22331/q-2021-07-13-501 (2021). - 38.
Kennes, D. M. et al. Moiréheterostructures as a condensed-matter quantum simulator.

*Nat. Phys.***17**, 155–163 (2021). - 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). - 40.
Georges, A., de’ Medici, L. & Mravlje, J. Strong correlations from hund’s coupling.

*Annu. Rev. Condens. Matter Phys.***4**, 137–178 (2013). - 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). - 42.
Mivehvar, F., Ritsch, H. & Piazza, F. Cavity-quantum-electrodynamical toolbox for quantum magnetism.

*Phys. Rev. Lett.***122**, 113603 (2019). - 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.
Melzi, R. et al. li

_{2}VO(*s**i*,*g**e*)*O*_{4}, a prototype of a two-dimensional frustrated quantum heisenberg antiferromagnet.*Phys. Rev. Lett.***85**, 1318–1321 (2000). - 45.
Kim, J. et al. Magnetic excitation spectra of sr

_{2}iro_{4}probed by resonant inelastic x-ray scattering: establishing links to cuprate superconductors.*Phys. Rev. Lett.***108**, 177003 (2012). - 46.
Banerjee, A. et al. Proximate kitaev quantum spin liquid behaviour in a honeycomb magnet.

*Nat. Mater.***15**, 733–740 (2016). - 47.
Jurcevic, P. et al. Quasiparticle engineering and entanglement propagation in a quantum many-body system.

*Nature***511**, 202–205 (2014). - 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). - 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). - 50.
Krešić, I. et al. Spontaneous light-mediated magnetism in cold atoms.

*Commun. Phys.***1**, 1–9 (2018). - 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). - 52.
Gao, H., Schlawin, F., Buzzi, M., Cavalleri, A. & Jaksch, D. Photoinduced electron pairing in a driven cavity.

*Phys. Rev. Lett.***125**, 053602 (2020). - 53.
Reschke, S., Mayr, F., Wang, Z., Do, Seung-Hwan, Choi, K.-Y. & Loidl, A. Electronic and phonon excitations in

*α*− RuCl_{3}.*Phys. Rev. B***96**, 165120 (2017). - 54.
Cao, H. B. et al. Low-temperature crystal and magnetic structure of

*α*− rucl_{3}.*Phys. Rev. B***93**, 134423 (2016).

## Acknowledgements

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.

## Funding

Open Access funding enabled and organized by Projekt DEAL.

## Author information

### Affiliations

### Contributions

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

## Ethics declarations

### Competing interests

The authors declare no competing interests.

## Additional information

**Peer review information** *Nature 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 http://creativecommons.org/licenses/by/4.0/.

## About this article

### Cite this article

Chiocchetta, A., Kiese, D., Zelle, C.P. *et al.* Cavity-induced quantum spin liquids.
*Nat Commun* **12, **5901 (2021). https://doi.org/10.1038/s41467-021-26076-3

Received:

Accepted:

Published:

DOI: https://doi.org/10.1038/s41467-021-26076-3

## Comments

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.