From pinch points to pinch lines: a new spin liquid on the pyrochlore lattice Owen Benton, 1 L. D. C. Jaubert, 1 Han Yan, 1 and Nic Shannon 1

1 Okinawa Institute of Science and Technology Graduate University, Onna-son, Okinawa 904-0395, Japan

The mathematics of gauge theories lies behind many of the most profound advances in physics in the last 200 years, from Maxwell's theory of electromagnetism to Einstein's theory of general relativity. More recently it has become clear that gauge theories also emerge in condensed matter, a prime example being the "spin ice" materials which host an emergent "electromagnetic" gauge field. In spin ice, the underlying gauge structure is revealed by the presence of "pinch-point" singularities in neutronscattering measurements. Here we report the discovery of a new spin liquid where the low-temperature physics is naturally described by the fluctuations of a tensor field with a continuous gauge freedom. This gauge structure underpins a novel form of spin correlations, giving rise to "pinch-line" singularities -line-like analogues of the pinch-points observed in spin ice. Remarkably, these features may already have been observed in the pyrochlore material Tb 2 Ti 2 O 7 .

Gauge symmetries are paramount in the understanding of many of the most fundamental theories of physics. Recent decades have seen an increasing appreciation of the role of gauge theories in condensed matter physics, emerging from the long-wavelength description of the collective behaviour of electrons. Emergent gauge theories have proved particularly important in the study of spin liquids -strongly fluctuating, disordered magnetic states, the description of which lies beyond the familiar territory of Landau theory [1][2][3][4][5].

The use of a gauge theory to describe the fluctuations of a spin liquid is exemplified by the case of the spin ice materials R 2 M 2 O 7 (R=Ho, Dy, M=Ti, Sn) [START_REF] Castelnovo | Magnetic monopoles in spin ice[END_REF][START_REF] Castelnovo | Spin ice, fractionalization and topological order[END_REF]. At low temperatures, the spin configurations in a spin ice are subject to a constraint directly analogous to Gauss' law for a magnetic field and consequently may be described in terms of a gauge theory. Among the many striking consequences of this is the observation of pinch-point singularities [Fig. 1(a), left panel] in the magnetic neutron scattering structure factor [START_REF] Henley | Power-law spin correlations in pyrochlore antiferromagnets[END_REF], as observed in Ho 2 Ti 2 O 7 [START_REF] Fennell | Magnetic Coulomb phase in the spin ice Ho2Ti2O7[END_REF]. Pinch-point scattering has also been observed in the putative quantum spin ice Tb 2 Ti 2 O 7 [START_REF] Fennell | Power-law spin correlations in the pyrochlore antiferromagnet Tb2Ti2O7[END_REF][START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF][START_REF] Guitteny | Anisotropic propagating excitations and quadrupolar effects in Tb2Ti2O7[END_REF]. However, in this case the experimental scattering shows pronounced butterfly-like features in the non spin-flip (NSF) channel and the scattering in the spin-flip (SF) channel shows narrow arm-like features extending along the h111i directions of reciprocal space, neither of which features are predicted for a spin ice. This raises the question of whether other types of spin liquid may be found amongst rare-earth pyrochlore magnets.

In this Article we introduce a new type of spin liquid on the pyrochlore lattice. This spin liquid arises on the phase diagram of a realistic model for pyrochlore magnets. As with spin ice, the theory of this spin liquid contains a gauge sym- metry. The nature of this theory is fundamentally different to the Maxwellian theory which describes spin ice, but just as the emergent gauge structure of spin ice reveals itself in pinch-point scattering, so the gauge structure of this new spin liquid has striking consequences for scattering experiments. We will show that at low temperatures, this gauge structure leads to line-like singularities along the h111i directions of reciprocal space, which we dub "pinch lines" since they are extended versions of the pinch points exhibited in spin ice. This is particularly interesting in the light of neutron scattering results on the pyrochlore magnets Tb 2 Ti 2 O 7 and Yb 2 Ti 2 O 7 which show strong, sharpening features along the h111i directions of reciprocal space. Indeed, our theory is able to account for several features of the diffuse scattering observed in Tb 2 Ti 2 O 7 [START_REF] Fennell | Power-law spin correlations in the pyrochlore antiferromagnet Tb2Ti2O7[END_REF][START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF][START_REF] Guitteny | Anisotropic propagating excitations and quadrupolar effects in Tb2Ti2O7[END_REF], which are unaccounted for by a theory based on a spin ice model -as shown in the right panels of Fig. 1.

Results

We begin with the most general, symmetry-allowed, Hamiltonian for nearest neighbour anisotropic exchange on the py- 1)], for parameters J1=J2=J4=0, J3 < 0. No phase transition is observed down to T = 0.001|J3|. Instead, the order parameter susceptibilities of neighbouring ordered phases exhibit a Curie-law crossover, characteristic of a Coulombic spin liquid [START_REF] Jaubert | Topological-sector fluctuations and Curie-law crossover in spin ice[END_REF]. The symbols used for different symmetry channels are shown in an inset.

rochlore lattice [START_REF] Curnoe | Quantum spin configurations in Tb2Ti2O7[END_REF][START_REF] Ross | Quantum excitations in quantum spin ice[END_REF][START_REF] Yan | Living on the edge: ground-state selection in quantum spin-ice pyrochlores[END_REF]:

H ex = X hiji S i • Jij • S j , J01 = 0 @ J 2 J 4 J 4 J 4 J 1 J 3 J 4 J 3 J 1 1 A (1)
where the exchange matrix J01 couples nearest neighbours along the r 01 = (0, 1, 1) direction and the other exchange matrices can be generated from J01 using point group operations. As shown in [START_REF] Yan | Living on the edge: ground-state selection in quantum spin-ice pyrochlores[END_REF][START_REF] Benton | Classical and quantum spin liquids on the pyrochlore lattice[END_REF], it is possible to map out the entire classical ground state phase diagram of Eq. ( 1) by an exact transcription of the Hamiltonian in terms of local fields defined on each tetrahedron [START_REF] Yan | Living on the edge: ground-state selection in quantum spin-ice pyrochlores[END_REF][START_REF] Benton | Classical and quantum spin liquids on the pyrochlore lattice[END_REF][START_REF] Mcclarty | Energetic selection of ordered states in a model of the Er2Ti2O7 frustrated pyrochlore XY antiferromagnet[END_REF]:

H ex = 1 2 X tet  A2 m 2 A2 + + E m 2 E + T2 m 2 T2 + T1 m 2 T1 + T1+ m 2 T1+ + E 0 (2)
where all the coefficients ↵ 0, E 0 is the ground state energy and the sum runs over all tetrahedra in the lattice.

The five fields m A2 , m E , m T2 , m T1 , m T1+ appearing in Eq. ( 2) are defined in the supplemental material. They transform according to the A 2 , E, T 2 , T 1 irreducible representations of the point group and have respective dimension 1, 2, 3, 3 and 3. Along a line of points in parameter space the three ordered phases which respectively maximise the fields m T1 , m T2 , m E become degenerate. This line includes the point J 1 = J 2 = J 4 = 0, J 3 < 0 [cf. Fig. 2(a)]. For parameter sets along this line of points we have A2 = E = T1 = 0, A2 , T1+ > 0 and the Hamiltonian is given by

H CL ex = 1 2 X tet [ A2 m 2 A2 + T1+ m 2 T1+ ] + E 0 . ( 3
)
In a classical ground state of Eq. ( 3) it must be the case that

m A2 = 0, m T1+ = 0 (4)
for every tetrahedron in the lattice. All of the results derived in this paper flow from the implementation of these constraints. These provide an exact description of the classical ground states along the line in parameter space where the three phases in Fig. 2(a) are degenerate. Observing the consequences of these constraints does not, however, require precise fine tuning of the Hamiltonian to Eq. ( 3). These constraints will also dominate the physics at finite temperatures for any choice of parameters where The system thus remains in a disordered but highly-correlated state down to low temperature.

We can understand the correlations of the spin liquid from Eq. ( 4). The demand that the fields m A2 and m T1+ vanish everywhere leaves the fields {m E , m T2 , m T1 } with freedom to fluctuate in the ground state. The spatial variation of these fluctuations is constrained by the fact that neighbouring tetrahedra share a spin, therefore a fluctuation of the local fields on one tetrahedron affects the values of the local fields on the neighbouring tetrahedra. The fields m E , m T2 , m T1 must therefore fluctuate in a correlated manner in order to avoid inducing violations of Eq. ( 4). In what follows we show how these correlated fluctuations can be understood in terms of the fluctuations of a tensor field with a continuous gauge freedom.

The constraints on the spatial variation of m E , m T2 , m T1 may be obtained from the continuity of fields between A and B sublattice tetrahedra. The ground state constraints [Eq. ( 4)] in fact imply a set of local conservation laws, on the lattice. These conservation laws in turn suggest that a coarse-graining approach can be successful in describing the fluctuations of m E , m T2 , m T1 . And, unlike the global conservation laws which underpin hydrodynamic theories, these fully local conservation laws can have consequences even for short wavelength fluctuations, as we shall see. Expanding the local constraints to leading order in a gradient expansion we find

r • m T1 = 0 (5) 0 B @ @ x m 1 E 1 2 @ y m 1 E + p 3 2 @ y m 2 E 1 2 @ z m 1 E p 3 2 @ z m 2 E 1 C A + p 3 2 r ⇥ m T2 3 2 sin( 0 T1 ) 0 @ @ y m z T1 + @ z m y T1 @ z m x T1 + @ x m z T1 @ x m y T1 + @ y m x T1 1 A = 0 (6)
where the angle 0 T1 is a function of the exchange parameters, defined in the supplemental materials. We wish to resolve the constraints ( 5)-( 6) naturally using a gauge-theoretic approach. Note that since m T1 appears in both constraints, we cannot simply introduce separate gauge fields to resolve Eqs. ( 5)- [START_REF] Castelnovo | Magnetic monopoles in spin ice[END_REF]. Instead, we incorporate the eight components of {m E , m T2 , m T1 } into a traceless tensor field B:

B = 0 B B @ m 1 E p 3 2 m z T2 3 sin( 0 T 1 ) 2 m z T1 p 3 2 m y T2 3 sin( 0 T 1 ) 2 m y T1 p 3 2 m z T2 3 sin( 0 T 1 ) 2 m z T1 1 2 m 1 E + p 3 2 m 2 E p 3 2 m x T2 3 sin( 0 T 1 ) 2 m x T1 p 3 2 m y T2 3 sin( 0 T 1 ) 2 m y T1 p 3 2 m x T2 3 sin( 0 T 1 ) 2 m x T1 1 2 m 1 E p 3 2 m 2 E 1 C C A (7)
Satisfaction of Eq. ( 6), along with the condition Tr[B] = 0 is guaranteed by the introduction of a symmetric, tensor field Y and writing

B = D • Y, D ⌘ 0 @ 0 @ z @ y @ z 0 @ x @ y @ x 0. 1 A (8)
The form of the matrix gauge field Y is then constrained by Eq. ( 5) which is satisfied if we take Y of the form

Y = 0 @ W z W y W z W x W y W x 1 A . (9)
We can generate alternative forms of Y by applying Abelian gauge transformations to Eq. ( 9) of the form

Y µ⌫ ! Y µ⌫ + @ µ @ ⌫ ⇣ (10)
The transformations of Eq. ( 10) leave the flux matrix B, and therefore the physical spin system, unchanged. The form of Y in Eq. ( 9) thus corresponds to a specific choice of gauge. The theory of our spin liquid is therefore invariant under a group of gauge transformations ⇣ 2 R. Abelian gauge transformations of a similar form to Eq. ( 10), acting on tensor fields also appear in the linearized theory of general relativity [19] and the theory of S = 2 gauge fields [START_REF] Fronsdal | Massless fields with integer spin[END_REF]. At low temperatures, where there are only fluctuations of the local fields m E , m T2 , m T1 , the free energy will be controlled by the entropy of these fluctuations. Coarse graining over some volume much larger than a unit cell but much smaller than the whole system, there will be more states available (and therefore more entropy) with small values of m E , m T2 , m T1 [4]. The most general symmetry-allowed Gaussian free energy describing small fluctuations of these fields, when written in terms of the tensor field Y, takes the form

F SL = T Z d 3 r V u.c.  a Tr[(D • Y) • (D • Y) T ] + b Tr[(D • Y) • (D • Y)] + c X µ (D • Y) 2 µµ ( 11
)
which is invariant under the gauge transformations of Eq. ( 10).

The distinctive nature of this spin liquid, and of the theory which describes it [Eq. [START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF]], can be revealed by neutron scattering experiments. This can be seen by calculating the correlation functions of the local fields m E , m T2 , m T1 in momentum space. In addition to displaying pinch point singularities at zone centres, these correlation functions are singular approaching any momentum q which is along the (h, h, h) directions of reciprocal space, or along any direction related by the lattice symmetry to (h, h, h). This contrasts with the case of the Coulombic spin liquid which occurs in the case of spin ice, where the correlation functions are only singular at the Brillouin zone centre. Since the fields m E , m T2 , m T1 are simply linear combinations of the spins, this singular behaviour will also show up in the spin structure factor S(q), measurable in neutron scattering experiments.

In the vicinity of one of these singularities, at T = 0, the scattering is approximated by

S(K + q k + q ? ) ⇡ X ↵ ↵ (K, q k ) q ↵ ? q ? q 2 ? ( 12
)
where K is a reciprocal lattice vector, q k is parallel to a h111i 11)]. Singular features are visible where these planes intersect h111i directions, as indicated by the black circles in each panel. These "pinch-line" singularities, Eq. ( 12), are characteristic of the gauge structure of the new spin liquid. (c), (f) spin structure factor calculated in finite temperature Monte Carlo simulation, in the same regions of reciprocal space. The pinch lines appear in the simulation results as sharp features around the point where h111i directions intersect the plane. (d), (g) a calculation of the structure factor made with a lattice based 1/N theory, described in the supplemental material, also exhibiting pinch-line singularities. The simulations were performed at a temperature T = 0.001 K for a cluster of N = 256000 sites and dimensions 40a0 ⇥ 40a0 ⇥ 10a0 where a0 is the linear size of a cubic unit cell. Results were calculated for parameters J1=0.042, J2=0.122, J3= 0.118, J4= 0.04 meV, with anisotropic g-tensor g ? /g k =1/3, in approximate correspondence to Tb2Ti2O7 [START_REF] Cao | Ising versus XY anisotropy in frustrated R2Ti2O7 compounds as "seen" by polarized neutrons[END_REF]. Since the crystal field ground state in Tb2Ti2O7 is a non-Kramers doublet, the finite value of g ? should be thought of as coming from mixing with the low lying crystal field excitation.

direction and q ? is orthogonal to that direction. The coefficients ↵ determine the "orientation" of the singularity in q-space. Their dependence on the Brillouin zone K may be thought of as a form factor determining the contribution of the fluctuations of each field m to the scattering in each Brillouin zone. The dependence on q k is smooth and near a zone center K one may write

↵ (K, q k ) ⇡ ↵ (K, 0) + O(q 2 k ).

For

↵ ⇢ / ↵ the structure factor in the limit q ? ! 0 will depend on the direction of approach and we have a singularity, along the entire h111i direction. Eq. ( 12) has the form of a pinch-point singularity extended into a line. We therefore will refer to it as a "pinch-line" singularity.

These "pinch lines" can be observed by taking planar cuts through scattering, which intersect these lines away from reciprocal lattice vectors [Fig. 3(a)]. This is illustrated using a T = 0 calculation of S(q) from the gauge theory [Eq. [START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF]] ] in Fig. 3(b), (e). For comparison, we show in Fig. 3(c), (f), the same quantity calculated at finite temperature within classical Monte Carlo simulation.

The simulation results show sharp features in the structure factor approaching the h111i directions, as predicted by the theory [Eq. [START_REF] Fennell | Power-law spin correlations in the pyrochlore antiferromagnet Tb2Ti2O7[END_REF]]. There is a small broadening of these singularities, coming from the finite temperature thermal fluctuations present in the Monte Carlo simulation. These features are even more clearly visible in the correlation functions of the local fields {m E , m T2 , m T1 }, shown in the supplemental materials. The presence of the pinch lines in the simulation results is a strong validation of our theory of the spin liquid regime.

The continuum theory [Eq. [START_REF] Fennell | Power-law spin correlations in the pyrochlore antiferromagnet Tb2Ti2O7[END_REF]] was derived from local constraints, with associated local conservation laws, and the structure of the theory is inherited from the structure of those local constraints. This leads us to expect that the pinch-line singularities will be robust features of the spin liquid, even at short wavelengths. We have confirmed this expectation using two independent, lattice-based calculations. Firstly, the sharpening of the scattering around the h111i directions is clearly seen in the Monte Carlo simulations in Fig. 3(c),(f).

Secondly, we have also performed a 1/N calculation of the spin correlations along the lines of that performed for the Heisenberg model in Ref. [START_REF] Isakov | Dipolar Spin Correlations in Classical Pyrochlore Magnets[END_REF]. This calculation also predicts pinch-line singularities along the h111i directions of reciprocal space, as shown in Fig. 3(d),(g). It is therefore apparent that these singularities are a robust feature of the spin liquid, arising from the structure of its ground-state constraints, which is captured by the continuum theory derived in this Article.

Discussion

Thus far we have uncovered a spin liquid described by a tensor field carrying a continuous gauge symmetry, arising in a particular limit of a realistic model for magnetism on the pyrochlore lattice [Eq. ( 1)]. The signal feature of this spin liquid is sharp line-like singularities along h111i directions of reciprocal space, which occur in addition to pinch point singularities at zone centers. These "pinch-line" singularities are unique to the spin liquid discussed in this paper and as such provide a very discriminating "smoking gun" signature of this novel magnetic state. In the light of this discovery it is interesting to consider two known pyrochlore materials, which are often discussed in the context of spin liquid physics: Tb 2 Ti 2 O 7 and Yb 2 Ti 2 O 7 .

Tb 2 Ti 2 O 7 has long been a focal point for discussion of 3-dimensional spin liquid physics [START_REF] Gardner | Cooperative paramagnetism in the geometrically frustrated pyrochlore antiferromagnet in Tb2Ti2O7[END_REF][START_REF] Molavian | Dynamically induced frustration as a route to a quantum spin ice state in Tb2Ti2O7 via virtual crystal field excitations and quantum many-body effects[END_REF][START_REF] Gardner | Magnetic pyrochlore oxides[END_REF]. While Eq. ( 1 As a point of comparison to these experiments, the behaviour of the structure factor S(q) in the spin flip (SF) and non spin flip (NSF) channels, appropriate to a polarised neutron scattering experiment with initial polarisation n||(1, 1, 0), is shown in Fig. 1, for the same set of exchange parameters as in Fig. 3. Narrow prominences are visible in the SF channel along the h111i directions [Fig. 1(a)]. There are also pinch points in both channels at Brillouin zone centers. The distribution, orientation and polarisation dependence of the pinch points observed in [START_REF] Fennell | Power-law spin correlations in the pyrochlore antiferromagnet Tb2Ti2O7[END_REF] is the same as that in Fig. 1. In particular, we are able to reproduce the shape of the features in the NSF channel, something which cannot be done with a spin ice based description. The possibility that the theory described in this Article could apply to Tb 2 Ti 2 O 7 is lent weight by a recent attempt at parameterizing a pseudo-spin Hamiltonian for Tb 2+x Ti 2 x O 7+y [START_REF] Takatsu | Quadrupole order in the frustrated pyrochlore Tb2+xTi2 xO7+y[END_REF] which places it close to the three-way phase boundary at which this spin liquid emerges in our classical treatment.

Spin liquid behaviour at finite temperature does not rule out the possibility of an magnetic order at lower temperature. Indeed, recent experiments have demonstrated the pres-ence of competing ordering phenomena in Tb 2 Ti 2 O 7 , with quadrupolar [START_REF] Taniguchi | Long-range order and spin-liquid states of polycrystalline Tb2+xTi2 xO7+y[END_REF][START_REF] Kadowaki | Composite Spin and Quadrupole Wave in the Ordered Phase of Tb2+xTi2 xO7+y[END_REF] and short range ordered antiferromagnetic states [START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF][START_REF] Kermarrec | Gapped and gapless short range ordered magnetic states with ( 1 2 , 1 2 , 1 2 ) wavectors in the pyrochlore magnet Tb2+xTi2 xO7+[END_REF][START_REF] Guitteny | Mesoscopic correlations in Tb2Ti2O7[END_REF] being observed depending on the sample stoichiometry and experimental cooling protocol. This is consistent with the nature of the spin liquid considered in our manuscript, which sits at the confluence of many competing orders. In particular, we note that the ground-state manifold of the spin liquid contains states consistent with the q ⇤ = (±1/2, ±1/2, ±1/2) order observed under field cooled conditions [START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF][START_REF] Kermarrec | Gapped and gapless short range ordered magnetic states with ( 1 2 , 1 2 , 1 2 ) wavectors in the pyrochlore magnet Tb2+xTi2 xO7+[END_REF]. These states can only be connected to the other states of the spin liquid by rotation of an O(L) number of spins, where L is the linear size of the system. This may suggest an explanation for the sensitivity to how the system is cooled -namely that field cooling may drive the system into a state from which it is hard to reach the other parts of the ground state manifold.

The combination of spin liquid physics and prominent features in the scattering along q k (1, 1, 1) is also strongly reminiscent of the discussion surrounding another pyrochlore: Yb 2 Ti 2 O 7 [START_REF] Bonville | Transitions and spin dynamics at very low temperature in the pyrochlores Yb2Ti2O7 and Gd2Sn2O7[END_REF][START_REF] Ross | Two-dimensional kagome correlations and field induced order in the ferromagnetic XY pyrochlore Yb2Ti2O7[END_REF][START_REF] Thompson | Rods of neutron scattering intensity in Yb2Ti2O7: compelling evidence for significant anisotropic exchange in a magnetic pyrochlore oxide[END_REF][START_REF] Chang | Higgs transition from a magnetic Coulomb liquid to a ferromagnet in Yb2Ti2O7[END_REF][START_REF] Yaouanc | Influence of short range correlations on the µSR polarization functions in the slow dynamic limit: application to the quantum spin-liquid system Yb2Ti2O7[END_REF]. Indeed, it has recently been argued that the unusual physics of this material springs from competition between the E and T 1 regions of the phase diagram in Fig. 2(a) [START_REF] Yan | Living on the edge: ground-state selection in quantum spin-ice pyrochlores[END_REF][START_REF] Jaubert | Are multiphase competition & order-by-disorder the keys to understanding Yb2Ti2O7?[END_REF]. In this context it is not unreasonable to imagine that the physics of the paramagnetic phase of Yb 2 Ti 2 O 7 may be influenced by a nearby spin liquid phase of the form described here. This provides an interesting alternative scenario to "quantum spin ice" physics in that material.

The theory presented in this Article provides a fundamentally different paradigm to the "emergent electromagnetism" known from spin ice and possesses a gauge freedom bearing an intriguing similarity to that appearing in the linearized theory of general relativity. This leads to the possibility of a unified theory of classical spin liquids on the pyrochlore lattice, and a classification of the above based on their associated gauge freedoms and the consequent singularities in their correlation functions. These issues will be explored further elsewhere.

The discovery of a new classical spin liquid is also a promising starting point for the discovery of new quantum spin liquid states [START_REF] Moessner | Three-dimensional resonatingvalence-bond liquids and their excitations[END_REF][START_REF] Hermele | Pyrochlore photons: The U(1) spin liquid in a S=1/2 three-dimensional frustrated magnet[END_REF]. This possibility is exemplified by the quantum spin ice problem, in which introducing quantum fluctuations amongst the manifold of classical spin ice ground states leads to a quantum spin liquid with dynamic U (1) gauge fields and fractionalized charges as excitations [START_REF] Hermele | Pyrochlore photons: The U(1) spin liquid in a S=1/2 three-dimensional frustrated magnet[END_REF]. Indeed, the limit of perturbative quantum fluctuations around a classical spin ice is currently the only limit in which a quantum spin liquid is reliably predicted by a variety of theoretical methods in a realistic model for a pyrochlore magnet [START_REF] Banerjee | Unusual liquid state of hard-core bosons on the pyrochlore lattice[END_REF][START_REF] Benton | Seeing the light: Experimental signatures of emergent electromagnetism in a quantum spin ice[END_REF][START_REF] Shannon | Quantum ice: A quantum Monte Carlo study[END_REF][START_REF] Savary | Coulombic quantum liquids in spin-1/2 pyrochlores[END_REF][START_REF] Hao | Bosonic many-body theory of quantum spin ice[END_REF][START_REF] Kato | Numerical evidence of quantum melting of spin ice: quantum-to-classical crossover[END_REF]. In the quantum spin ice case the classical and quantum limits of the problem are linked by a smooth thermal crossover [START_REF] Benton | Seeing the light: Experimental signatures of emergent electromagnetism in a quantum spin ice[END_REF][START_REF] Kato | Numerical evidence of quantum melting of spin ice: quantum-to-classical crossover[END_REF]. Similarly, studies of the quantum S = 1/2 Heisenberg model on the pyrochlore lattice find similar spin correlations [START_REF] Canals | Pyrochlore Antiferromagnet: A Three-Dimensional Quantum Spin Liquid[END_REF][START_REF] Canals | Quantum spin liquid: The Heisenberg antiferromagnet on the three-dimensional pyrochlore lattice[END_REF][START_REF] Huang | Spin-Ice State of the Quantum Heisenberg Antiferromagnet on the Pyrochlore Lattice[END_REF] to those predicted by the gauge theoretic description of the classical problem [START_REF] Isakov | Dipolar Spin Correlations in Classical Pyrochlore Magnets[END_REF]. We therefore expect that the classical results presented for the spin liquid discussed here are likely to remain relevant upon the inclusion of quantum fluctuations.

In conclusion, we have demonstrated the existence of a new classical spin liquid phase on the pyrochlore lattice, described by the fluctuations of a tensor field with a continuous gauge freedom. The unusual gauge structure is revealed by "pinchline" singularities in correlations which could be observed in neutron scattering experiments.

From pinch points to pinch lines: a new spin liquid on the pyrochlore lattice-supplemental material Owen Benton, 1 Ludovic Jaubert, 1 Han Yan, 1 and Nic Shannon 1

1 Okinawa Institute of Science and Technology Graduate University, Onna-son, Okinawa 904-0395, Japan

DEFINITION OF LOCAL FIELDS AND SUSCEPTIBILITIES

It was shown in [1] that the generalized model for nearest neighbour exchange on the pyrochlore lattice [Eq. (1) of the main text] may be exactly rewritten in terms of local fields, defined on the pyrchlore tetrahedra

H ex = 1 2 X tet A2 m 2 A2 + E m 2 E + T2 m 2 T2 + T1+ m 2 T1+ + T1 m 2 T1 + constant. (1)
The fields are labelled by the irreducible representations of the point group according to which they transform. These fields are defined in Table I. (1) of the main text] reduces to a sum of quadratic forms when written in terms of these fields [1]. The labels A2, E, T2, T1 refer to the irreducible representation of the point group according to which the corresponding field transforms. There are two fields which transform according to the T1 irrep. The parameter 0 T 1 is chosen as a function of the exchange parameters to remove the symmetry allowed coupling between these fields.

Local OP Definition in terms of spins m

A 2 1 p 3 (S x 0 + S y 0 + S z 0 + S x 1 S y 1 S z 1 S x 2 + S y 2 S z 2 S x 3 S y 3 + S z 3 ) m E 1 p 6 ( 2S x 0 + S y 0 + S z 0 2S x 1 S y 1 S z 1 + 2S x 2 + S y 2 S z 2 + 2S x 3 S y 3 + S z 3 ) 1 p 2 ( S y 0 + S z 0 + S y 1 S z 1 S y 2 S z 2 + S y 3 + S z 3 ) ! m T 2 0 B B @ 1 p 2 ( S y 0 + S z 0 + S y 1 S z 1 + S y 2 + S z 2 S y 3 S z 3 ) 1 p 2 (S x 0 S z 0 S x 1 S z 1 S x 2 + S z 2 + S x 3 + S z 3 ) 1 p 2 ( S x 0 + S y 0 + S x 1 + S y 1 S x 2 S y 2 + S x 3 S y 3 ) 1 C C A m T 1,ice 0 B B @ 1 p 3 (S x 0 + S y 0 + S z 0 + S x 1 S y 1 S z 1 + S x 2 S y 2 + S z 2 + S x 3 + S y 3 S z 3 ) 1 p 3 (S x 0 + S y 0 + S z 0 S x 1 + S y 1 + S z 1 S x 2 + S y 2 S z 2 + S x 3 + S y 3 S z 3 ) 1 p 3 (S x 0 + S y 0 + S z 0 S x 1 + S y 1 + S z 1 + S x 2 S y 2 + S z 2 S x 3 S y 3 + S z 3 ) 1 C C A m T 1,planar 0 B B @ 1 p 6 ( 2S x 0 + S y 0 + S z 0 2S x 1 S y 1 S z 1 2S x 2 S y 2 + S z 2 2S x 3 + S y 3 S z 3 ) 1 p 6 (S x 0 2S y 0 + S z 0 S x 1 2S y 1 + S z 1 S x 2 2S y 2 S z 2 + S x 3 2S y 3 S z 3 ) 1 p 6 (S x 0 + S y 0 2S z 0 S x 1 + S y 1 2S z 1 + S x 2 S y 2 2S z 2 S x 3 S y 3 2S z 3 ) 1 C C A m T 1 cos( 0 T
The angle 0 T1 which appears in the definitions of m T1+ and m T1 is chosen such that there is no bilinear coupling between m T1+ and m T1 and such that a T1  a T1+ .

(

) 2
Note that this convention for the definition of the T 1 symmetric fields is dfifferent to that chosen in Ref. [1]. Following common practice in Monte Carlo simulations, the order parameter susceptibilities appearing in Fig. 2 of the main text are calculated in according to the following formula:

= N T hm 2 i h|m |i 2 (3)
where N is the number of spins in the system, T is the temperature and m are the local fields given in Table I.

(a) hm 1 E (q)m 1 E ( q)i (b) hm 1 E (q)m 1 E ( q)i (c) hm 1 E (q)m 1 E ( q)i (d) hm z

CORRELATION FUNCTIONS OF THE FIELDS m

In Figs. 1 we show examples of the momentum space correlation functions of the local fields {m E , m T2 , m T1 }, as calculated from the continuum theory developed in the manuscript (left half of panels) and from Monte Carlo simulation (right half of panels). Specifically, we show the correlation functions hm E (q) 1 m E ( q) 1 i, hm T2 (q) z m T2 ( q) z i, in three planes of reciprocal space q = 2⇡ a 0 (h, k, 0), q = 2⇡ a 0 (h, k, 0.2), q = 2⇡ a 0 (h, k, 0.4)

in each case we see singular features wherever the planes cut the h111i directions. The Fourier transforms of the local fields are defined as

m ↵ (q) = r 1 N u.c. X tet2A exp[ iq • r tet ]m ↵ (r tet ) (4)

FIG. 1 :

 1 FIG. 1: Comparison of correlations in spin ice (left half of panels) with those of the new spin liquid (right half of panels). Predictions for polarised neutron scattering experiments are shown in the (a) spin-flip (SF) and (b) non spin-flip (NSF) channels, as measured by Fennell et al. [9, 10]. The prediction for spin ice exhibits pinchpoint singularities in the SF channel. The NSF channel is completely featureless in a nearest neighbour model for spin ice-as shown hereand develops smooth maxima at the zone boundaries in the presence of long range dipole interactions [9]. In contrast, the new spin liquid exhibits singular features in both SF and NSF scattering. Results are taken from classical Monte Carlo simulation of the nearest neighbour model Hex [Eq. (1)], as described in the text.

1 FIG. 2 :

 12 FIG. 2: Evidence of spin-liquid behaviour from Monte Carlo simulation. (a) Classical ground state phase diagram of Hex [Eq. (1)] for J3<0, in the plane J4=0, showing how ordered phases with symmetry T1, E and T2 meet at the point J1=J2=0 [15]. (b) Order-parameter susceptibilities and (c) heat capacity calculated in classical Monte Carlo simulation of Hex [Eq. (1)], for parameters J1=J2=J4=0, J3 < 0. No phase transition is observed down to T = 0.001|J3|. Instead, the order parameter susceptibilities of neighbouring ordered phases exhibit a Curie-law crossover, characteristic of a Coulombic spin liquid[START_REF] Jaubert | Topological-sector fluctuations and Curie-law crossover in spin ice[END_REF]. The symbols used for different symmetry channels are shown in an inset.

  cost of having a finite value of the fields m T1 , m T2 , m E is much lower than the cost to have a finite value of m A2 , m T1+ . The constraints in Eq. 4 are insufficient to select an ordered ground state in themselves. In such circumstances, fluctuations may select a preferred ordered state via the order-bydisorder mechanism, but Monte Carlo simulations indicate that they fail to do so, down to temperatures 3 orders of magnitude below the scale of the bare coupling [see Fig. 2(b)].

FIG. 3 :

 3 FIG.3: Gauge structure of the new spin liquid, as revealed through "pinch-line" singularities. (a) Location of pinch-line singularities in reciprocal space. (b), (e) spin structure-factor S(q) in parallel planes in reciprocal space in the Brillouin zone centred on K = (3, 3, 1), as calculated from the continuum theory [Eq.[START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF]]. Singular features are visible where these planes intersect h111i directions, as indicated by the black circles in each panel. These "pinch-line" singularities, Eq. (12), are characteristic of the gauge structure of the new spin liquid. (c), (f) spin structure factor calculated in finite temperature Monte Carlo simulation, in the same regions of reciprocal space. The pinch lines appear in the simulation results as sharp features around the point where h111i directions intersect the plane. (d), (g) a calculation of the structure factor made with a lattice based 1/N theory, described in the supplemental material, also exhibiting pinch-line singularities. The simulations were performed at a temperature T = 0.001 K for a cluster of N = 256000 sites and dimensions 40a0 ⇥ 40a0 ⇥ 10a0 where a0 is the linear size of a cubic unit cell. Results were calculated for parameters J1=0.042, J2=0.122, J3= 0.118, J4= 0.04 meV, with anisotropic g-tensor g ? /g k =1/3, in approximate correspondence to Tb2Ti2O7[START_REF] Cao | Ising versus XY anisotropy in frustrated R2Ti2O7 compounds as "seen" by polarized neutrons[END_REF]. Since the crystal field ground state in Tb2Ti2O7 is a non-Kramers doublet, the finite value of g ? should be thought of as coming from mixing with the low lying crystal field excitation.

  ) alone may not constitute a complete quantitative model for the physics of Tb 2 Ti 2 O 7 it is interesting to compare observations on Tb 2 Ti 2 O 7 with the phenomenology of the spin liquid. Polarized neutron scattering experiments on Tb 2 Ti 2 O 7 have shown evidence of singular scattering at Brilllouin zone centers, but the form of this scattering looks rather different to a typical spin ice, especially in the non-spin flip (NSF) channel. At the same time, the data presented in Ref. 10, shows bright, narrow features extending along the h111i directions.

  Definition of the local fields m appearing in Eq. (1). The exchange Hamiltonian Hex [Eq.

FIG. 1 :

 1 FIG.1: Comparison of the correlation functions hm 1 E (q)m 1 E ( q)i, hm z T 2 (q)m z T 2 ( q)i as calculated in the continuum field theory (left half of panels) and Monte Carlo simulation (right half of panels) of a cluster of N = 256000 spins at temperature T = 0.001 K. The field theory calculation is shown on a discrete grid in momentum space for the sake of comparison with the simulations which are carried out on a finite lattice with linear dimensions Lx ⇥ Ly ⇥ Lz = 40 ⇥ 40 ⇥ 10 a 3 0 with periodic boundary conditions, where a0 is the linear dimension of a cubic unit cell. The singular features found in simulation are exactly reproduced by the continuum theory [Eq.[START_REF] Fritsch | Antiferromagnetic correlations at (1/2, 1/2, 1/2) in the ground state of the pyrochlore magnet Tb2Ti2O7[END_REF] of main text], with pinch-line singularities visible in both calculations where the planes of scattering cut the h111i directions. Calculations are carried out for the set of exchange parameters J1=0.042, J2=0.122, J3= 0.118, J4= 0.04 meV, as in Fig.3of the mauscript.

Acknowledgements

This work was supported by the Theory of Quantum Matter Unit of the Okinawa Institute of Science and Technology Graduate University.

Author contributions

The analytic calculations were performed by O. B. and H. Y. The Monte Carlo simulations were written and performed by L. D. C. J. The manuscript was written by O. B. with input from all authors. N. S. supervised the project.

where N u.c. is the number of unit cells, the sum runs over only the 'A' sublattice of tetrahedra and r tet are the positions of the centers of those tetrahedra.

The correlation functions agree well between the continuum theory and simulation. The qualitative structure of the correlations is in agreement for all cases, although some quantitative differences are visible at shorter wavelengths. Most importantly, the Monte Carlo simulations clearly show the sharpening of the correlation functions approaching the h111i directions, as predicted by the field theory, demonstrating that the pinch lines are a robust feature of the spin liquid regime.

LATTICE BASED CALCULATION OF THE CORRELATIONS IN THE SPIN-LIQUID REGIME

For the purposes of comparison with the continuum, we have also performed some fully lattice based calculations of the correlations in the spin liquid regime.

These calculations follow a method which has been previously been shown successful in understanding the correlations of disordered phases of the Heisenberg model on the pyrochlore lattice [2], spin ice [3] and protons in water ice [4].

In this approach the constraints on the lengths of the spins

are only enforced on average

Eq. ( 6) is enforced by means of a Lagrange multiplier added to the Hamiltonian. We write

where is the inverse temperature.

Using a Fourier transformation H may be written as

where is S(q) is a 12-component vector formed from the Fourier transforms of the 3 spin components on each of the 4 sublattices. The correlations of S(q) are then

and can be chosen such that Eq. ( 6) is obeyed.

Where M(q) possesses flat bands of eigenvalues at the bottom of its spectrum-as is the case in the spin liquid regimelimit T ! 0 of the correlation function becomes a projection matrix, projecting into the subspace described by the associated eigenvectors [4]. This projection operator can be thought of as enforcing the local ground state constraints [3].

It is this, zero temperature limit of the correlation function which is plotted in Figs. 3(d) and 3(g) of the main text. The approach outlined here can be constructed as a pertubative expansion in powers of 1/N , where N is a number of copies of the system and the spin length constraint [Eq. ( 6)] becomes

This method is described in more detail in Ref. [5].


L D C Jaubert

Monopole holes in a partially ordered spin liquid

If spin liquids have been famously defined by what they are not, i.e. ordered, the past years have seen the frontier between order and spin liquid starting to fade, with a growing number of materials whose low-temperature physics cannot be explained without co-existence of (partial) magnetic order and spin fluctuations. Here we study an example of such co-existence in the presence of magnetic dipolar interactions, related to spin ice, where the order is long range and the fluctuations support a Coulomb gauge field. Topological defects are effectively coupled via energetic and entropic Coulomb interactions, the latter one being stronger than for the spin-ice ground state. Depending on whether these defects break the divergence-free condition of the Coulomb gauge field or the long-range order, they are respectively categorized as monopoles -as in spin ice -or monopole holes, in analogy with electron holes in semiconductors. The long-range order plays the role of a fully-occupied valence band, while the Coulomb spin liquid can be seen as an empty conducting band. These results are discussed in the context of other lattices and models which support a similar co-existence of Coulomb gauge field and long-range order. We conclude this work by explaining how dipolar interactions lift the spin liquid degeneracy at very low energy scale by maximizing the number of flippable plaquettes, in light of the equivalent quantum dimer model.

The possibility to recast the collective behavior of electrons and atoms as elegant emergent phenomena is probably one of the most fascinating aspect of condensed matter. The emergence of quasi-particles does not only allow for a deeper understanding of the problem at hand but has also often led to surprising connections across physics, as recently illustrated by photon-like magnetic excitations in quantum spin ice [1][2][3][4] . Such approach takes an enhanced flavor when the particle has not yet been observed at high energy, or may not exist. This is for example the case for Majorana fermions observed in nanowires coupled to superconductors 5 , and for magnetic monopoles and their underlying Coulomb gauge field 6 .

Coulomb gauge theories have emerged from a variety of discrete models, on the kagome lattice 7 , in fully-packed loop models 8 , itinerant-electron systems at partial-filling [9][10][11] and chemically disordered materials such as CsNiCrF 6 12 . We refer the interested reader to the excellent review by Chris Henley on this topic 13 . But its most famous experimental realization has probably been observed in the classical spin liquid ground state of spin ice materials Dy 2 Ti 2 O 7 and Ho 2 Ti 2 O 7 , where topological excitations take the form of magnetic monopoles effectively interacting via long-range Coulomb interactions.

Interestingly, a Coulomb phase is not incompatible with partial magnetic order. For example, the stability of a Coulomb ferromagnet has been discussed in the context of the generalized Quantum Spin Ice model 14,15 and observed over a finite temperature window in a classical spin ice model with broken rotational symmetry 16 . In two dimensions, such co-existence has been studied in the two-stage ordering process of the dipolar kagome ice model 17,18 , of direct interest for artificial magnetic lattices 19 . In higher dimensions, advanced numerical simulations of the nearest-neighbour valence-bond wave functions on the cubic and diamond lattices have shown persistence of a small but clearly finite squared staggered moment in the thermodynamic limit 20 .

Our goal in this paper is to study what happens when a Coulomb spin liquid co-exists with long-range order on the pyrochlore lattice [21][22][23] . In particular, based on the Helmholtz decomposition of Fig. 2, we explain in section V how quasi-particles which break the long-range order can be understood as monopole holes in analogy with electron holes in semiconductors. In that sense, the Coulomb spin liquid plays the role of a conduction band where defects are "standard" monopoles, while the long-range ordered part of the spin degrees of freedom serves as a fully occupied valence band.

After a brief presentation of the physics of spin ice in section I, we introduce the concept of a Fragmented Coulomb Spin Liquid (FCSL) in section II which will be the central theme of our paper. In particular, the topology of the FCSL and its mapping onto the hard-core dimer model on the diamond lattice are discussed, before presenting the numerical methods used in this work. The nearest-neighbour version of the FCSL is briefly studied in section III in order to compute the strength of the entropic Coulomb potential between topological defects. Magnetic dipolar interactions are included in section IV and the FCSL is explained in terms of the dumbbell model 6 , focusing on the properties due to the underlying co-existence with long-range order. Notably, topological defects can be divided into two kinds of quasiparticles and are found to interact via an effective energetic Coulomb potential. Based on these results, section V is devoted to the similitudes between the FCSL and semiconductor physics. The relevance of our work to other systems, especially on the cubic lattice, is discussed in section V B. To conclude, in section VI, we confirm that dipolar interactions order the FCSL at very low temperature 21 whose ground state we identify to be the so-called R-states that have been studied in the quantum dimer model on the diamond lattice [24][25][26] . Additionally, the pinch points in the structure factor of Fig. 10 are found to persist down to the transition temperature, despite the build-up of R-states correlations.

I. WHAT IS THE COULOMB SPIN LIQUID IN SPIN ICE ?

The spin ice model consists of Ising spins S i on the pyrochlore lattice made of corner-sharing tetrahedra. The pyrochlore structure is given in Fig. 1. Each spin belongs to two tetrahedra and can only point towards the center of these two tetrahedra. In spin-ice materials, this constraint is imposed by the surrounding crystal field of the oxygens. Thus for a given tetrahedron, its four spins are not collinear and form four distinct sublattices. For nearest-neighbour ferromagnetic interactions, the ground state of this model is a highly degenerate classical spin liquid 27,28 where all tetrahedra have two spins pointing in and two spins pointing out. The resulting coarse-grained magnetization field is locally conserved and divergencefree, conferring the name of Coulomb spin liquid to the ground state of spin ice, in analogy to Maxwell's equations. Experimentally, the resulting spin-spin correlations are dipolar-like and take the form of pinch-points in the structure factor, as measured by neutron scattering in Ho 2 Ti 2 O 7 [29][30][31][32] . But what happens beyond nearest-neighbours ? The question is especially relevant to magnetic dipolar interactions (see Eq. ( 1)), which can be important for rareearth ions where nearest-neighbour couplings are rather weak (∼ 1K) and single-ion magnetic moments very large (up to ∼ 10µ B ):

E dip = Dr 3 p i>j S i • S j -3 S i • e ij S j • e ij r 3 ij ( 1
)
where i, j are pyrochlore sites separated by a distance r ij along the unit-length vector e ij and r p is the nearestneighbour distance on the pyrochlore lattice. The size of the magnetic moment µ is included in the prefactor

D = µ 0 µ 2 4π r 3 p (2)
where µ 0 is the vacuum magnetic permeability. Remarkably, the Coulomb-spin-liquid degeneracy is only weakly lifted by such long-range interactions 33,34 , because magnetic dipolar interactions between spins can be recast as effective magnetic Coulomb interactions between topological excitations, namely the magnetic monopoles 6,35,36 . The spin-ice ground state thus appears as a vacuum of magnetic charges. Following the procedure developed in Ref. [6], the Coulomb potential between two magnetic charges at distance r is

V = - µ 0 Q 2 4π r where Q = ± 2µ r d (3)
where r d is the nearest-neighbour distance on the diamond lattice and shortest distance between charges. Using Eq. ( 2), the potential energy can be expressed dimensionless as follow

V D = - 8 3 2 3 1 r/r d (4)
At very low energy scales, the quasi-degeneracy of the Coulomb spin liquid is ultimately lifted by dipolar interactions [37][38][39] and leads to exotic physics in a magnetic field, such as the magnetisation plateau in a [001] field analogue to a quantum solid phase in (2 + 1) dimensions 40 .

II. WHAT IS A FRAGMENTED COULOMB SPIN LIQUID ?

A. Presentation

In spin ice, the nearest-neighbour coupling plays the role of an effective chemical potential for monopoles 6 . Using chemical pressure 42,43 , it is possible to tune the FIG. 1. The Fragmented Coulomb Spin Liquid (FCSL) studied here is made of Ising spins on the pyrochlore lattice where all up (blue) and down (red) tetrahedra respectively carry effective positive ("3 in -1 out") and negative ("3 out -1 in") magnetic charges. The magnetic charges crystallize in a zinc-blende structure on the diamond lattice, whose sites are in the centers of the tetrahedra. Time-reversal symmetry reverses the sign of all charges, but conserves the FCSL ensemble. The minority spins are colored in green and shared by both neighbouring tetrahedra. A given spin configuration is entirely determined by the minority spins, modulo timereversal symmetry, which provides a two-to-one mapping onto a hard-core dimer model on the diamond lattice 24,26,41 . Please note that the fully saturated configuration in a [111] magnetic field is part of the FCSL ensemble. FIG. 2. The Helmholtz decomposition of the spin degrees of freedom is done as follow 22 . Each spin contributes to two magnetization fluxes. Left: for the double charges "4 in" (or "4 out"), all fluxes are pointing inwards (or outwards), giving rise to maximum or full divergence. Right: for the "2 in -2 out" configurations, there are as many fluxes pointing inwards as outwards, giving rise to a divergence-free field, i.e. with local flux conservation. Center : for the "3 in -1 out" configurations, the two contributions are noticeably different, being the sum of a full -and a zero -divergence. In the FCSL, the former contribution represents the long-range charge order, while the latter contribution is the Coulomb spin liquid. Creating a "4 in / 4 out" or a "2 in -2 out" defect out of the FCSL locally breaks respectively the zeroor the full-divergence. chemical potential and favor monopoles. Theoretically, if we forbid the presence of "4 in / 4 out" tetrahedra, a large negative chemical potential shall order the singlycharged monopoles into a zinc-blende structure 21,22 , as depicted on Fig. 1. In this long-range charge order, all positive charges have four negatively charged nearest neighbours, and vice-versa. As illustrated in Fig. 2, the Helmholtz decomposition of the spin degrees of freedom 22 shows that a fraction of the degrees of freedom forms a divergence-full configuration (the zinc-blende charge order) while the remaining part is divergence-free and supports an extensive degeneracy. In this phase, the local divergence condition is thus half-full, or arguably, half-empty. To emphasize this feature, we shall refer to this partially ordered phase as a Fragmented Coulomb Spin Liquid (FCSL). In particular, we use the term "fragmented" as a generalization since other phases may present the same physics with various degree of magnetic order (see discussion in section V B).

The competition between the FCSL and the spinice ground state has been studied for fixed density of monopoles 21 and fixed chemical potential 22 in absence of "4 in / 4 out" tetrahedra. However, in order to study the general properties of the Coulomb spin liquid coexisting with magnetic order, we do not want to eliminate any configuration a priori. In presence of "4 in / 4 out" tetrahedra, the FCSL is a priori unstable with only nearest-neighbour Ising spin-spin coupling and dipolar interactions 22,23 . But it can be stabilized for example via four-body interactions

H = J ν   ij ∈ν S i • S j   2 , (5)
where ν runs over all tetrahedra in the system and ij represents the six pairs of spins in tetrahedron ν. As illustrated in Fig. 3 for Hamiltonian (5), the eight single charges are in the ground state ("3 in -1 out" and "3 out -1 in"), while the six "2 in -2 out" and two "4 in / 4 out" configurations are excited states. Dipolar interactions would then naturally order the singly-charged monopoles into a zinc-blende structure as in Fig. 1. Further nearest-neighbour interactions can then easily tune the energy cost of the "2 in -2 out" and "4 in / 4 out" configurations. It should also be noted that the FCSL can be partially stabilized for a precise set of parameters of the dipolar spin ice model 23 . More generally, singlycharged monopoles may be favored by quantum fluctuations in the context of Yb 2 Ti 2 O 7 44 and magneto-electric coupling 45,46 .

However, our goal in this paper is not to focus on a given Hamiltonian or to propose a mechanism to stabilize the FCSL phase, which will be investigated elsewhere 47 . Our goal is to describe and understand how the co-existence of order and spin liquid, mediated by dipolar interactions, affect the emergence of topological defects and eventually the ordering of the spin liquid phase itself. From now on, we shall assume that the FCSL is stable and that "2 in -2 out" and "4 in / 4 out" spin configurations are excitations out of the FCSL.

B. Pseudo-magnetization

A convenient way to distinguish between the different phases that will be encountered in this work is via the pseudo-magnetization

ρ = 1 N i S i • e i = 1 N i σ i (6)
where e i is the local easy-axis vector pointing out of down tetrahedra. Hence the pseudo-spin σ i = S i • e i = -1 if the spin points inwards a down tetrahedron, and +1 if it points outwards. The "2 in -2 out" spin-ice ground state, FCSL and antiferromagnetic "all in / all out" long-range order respectively take the values of ρ = {0, 1/2, 1}, both globally for the entire system and locally on each tetrahedron. This quantity will be of particular interest when discussing the semiconductor analogy in section V.

FIG. 3. The fragmented Coulomb spin liquid considered in this paper is entirely composed of single charges ("3 in -1 out" and "3 out -1 in"), as opposed to the spin-ice ground state made of "2 in -2 out" tetrahedra. The "2 in -2 out" and "4 in / 4 out" configurations are considered as excitations out of the FCSL. All blue (resp. red) spins are pointing inwards (resp. outwards) of these up tetrahedra.

C. The dimer model on the diamond lattice

The minority spins of the FCSL can be recast as dimers on the diamond lattice, as shown in Fig. 1 where the dimers are the green spins. This mapping allows for a direct estimate of the degeneracy of the FCSL, with approximately 1.3 N/2 configurations 41 , where N is the total number of pyrochlore sites. This value is to be compared to the spin-ice ground state degeneracy, with approximately 1.5 N/2 configurations 48,49 . The nearestneighbour version of our model is thus the classical version of the quantum dimer model on the diamond lattice [24][25][26] , which has been discussed in the context of the magnetization plateau observed in HgCr 2 O 4 and CdCr 2 O 4 [50][51][52] . This parallel will be useful in section VI.

D. Topology

Starting from a given FCSL configuration and flipping a closed chain made of spins which are alternatively pointing in and out of successive neighboring tetrahedra, automatically creates another FCSL configuration. This is because such a spin update respects both the full-and zero-divergence constraints of the FCSL phase. This closed chain of spins is often better known as a "worm" because it is used as an update of the worm algorithm in Monte Carlo simulations. A worm can be local and as small as 6 sites, as illustrated by the thick lines on Fig. 9, or on the contrary extensively long and winding around the system with periodic boundary conditions. The "2 in -2 out" and "4 in / 4 out" excitations mentioned previously are actually topological defects whose creations, diffusion and annihilations correspond to worm updates.

The topology of the FCSL phase is thus strongly reminiscent of the topology of the "2 in -2 out" spin-ice ground state 13 , but with one main difference. In a topological phase, two configurations belonging to the same topological sector are connected by non-winding-worm updates. When winding worms are necessary and sufficient to transform one configuration into the other, the two configurations belong to different topological sectors, but to the same Kempe sector 53,54 . A trivial example is the fully saturated configuration along a given [111] direction, which is a topological sector by itself. When two configurations cannot be related by any worm update, then they are said to belong to different Kempe sectors. In the "2 in -2 out" ensemble, any pair of configurations can be transformed into each other by successive worm updates; there is thus only one Kempe sector. This is not the case for the FCSL. By construction, a worm cannot transform a positive charge into a negative one, or viceversa. Thus two FCSL configurations connected by timereversal symmetry cannot be connected by worms. It is easy to show that the FCSL is actually divided into two different Kempe sectors. As a consequence of the timereversal symmetry, a given Kempe sector is selected by a saturating magnetic field in the [111] direction, while reversing the magnetic field selects the other Kempe sector. It means there is a bijective mapping between the dimer model on the diamond lattice and each of the FCSL Kempe sectors. In absence of any time-reversal symmetry breaking term, the properties of the two Kempe sectors are the same.

E. Numerical methods

The dipolar energy of Eq. ( 1) has been implemented by Ewald summation into classical Monte Carlo simulations, without a demagnetization factor 55,56 . A variation of the worm algorithm 22,56 has been used to separate topological defects in sections III and IV, and to equilibrate the system in section VI. In the latter case, acceptance of the worm updates was done via a Metropolis argument. To further help thermalization in section VI, parallel tempering was also included 57,[START_REF] Geyer | Computing Science and Statistics: Proceedings of the 23rd Symposium on the Interface[END_REF] .

III. ENTROPIC COULOMB INTERACTIONS IN THE NEAREST-NEIGHBOUR MODEL

Even though one of the foci of the present paper is on the influence of long-range dipolar interactions, we shall briefly consider nearest-neighbour physics in this section.

Due to the underlying divergence-free condition, Coulomb phases are known to support entropic dipolar correlations between spins 13,[START_REF] Isakov | [END_REF]60 , which can be recast as an entropic interaction between topological defects. In the case of the "2 in -2 out" spin-ice phase, this entropic interaction has been computed analytically and checked numerically by Castelnovo et al. 61 , and found to be

E ent = α SI k B T R
, where α SI = 0.36755 (7) where T is the temperature and R = r/r d is the dimensionless distance between topological defects in units of the diamond nearest-neighbour distance r d . As a consequence, the probability for two topological defects in a Coulomb phase to be separated by a vector R is independent of the temperature

P ( R) ∝ exp E ent /k B T = exp α/| R| . ( 8
)
To measure such distribution is not only the signature of an emergent Coulomb gauge theory but also a way to characterize the Coulomb phase by the value of its parameter α. Since the exponential law of Eq. ( 8) is independent of the temperature, it remains valid in the limit T → 0 + where we can assume that only one pair of topological defects persists. This is why we have performed Monte Carlo simulations to measure the probability distribution for a unique pair of topological defects to be separated by vector R, or equivalently, the probability distribution to be separated by distance R = | R|, renormalized by the number of diamond sites sitting on a sphere of radius R. Technically, starting from a random FCSL configuration, we create a pair of topological defects which we separate using a worm algorithm. At each step, one of the defect is randomly moved to a neighbouring diamond site, keeping the number of defects constant, until the pair of defects is annihilated. The distribution is updated at each step. The process is repeated over 10 4 closed worms and averaged over 100 independent initial configurations. In order to check the validity of the exponential law over large distances, we considered a system made of 16 millions (1.6 10 7 ) sites. To validate our simulations against existing literature, we have also simulated the "2 in -2 out" spin-ice phase.

As shown in Fig. 4, our numerical results for spin ice (green triangles) perfectly match the theoretical value of α SI = 0.36755 (solid black line). As for the FCSL, we find a higher entropic interaction of α FCSL = 0.473 ± 0.005. For the sake of completeness, we have considered separately the cases where the initial pair of defects out of the FCSL were "2 in -2 out" and "4 in / 4 out", but the parameter α FCSL was found to be the same in both cases. The stiffness of a Coulomb phase is inversely related to the variance of the divergence-free field 60 and is known to be proportional to α 61 . Since the coarse-grained magnetization field of the FCSL is not divergence-free, one needs to be careful when defining the stiffness of the FCSL. However, since the divergence-full field is antiferromagnetic, it does contribute to the magnetization of the system. Magnetic fluctuations of the FCSL are FIG. 4. Monte Carlo simulations confirm that, in absence of dipolar interactions, the entropic probability distribution for a unique pair of topological defects to be separated by a vector R follows an exponential scaling law as in Eq. ( 8). This scaling law for the "2 in -2 out" spin-ice phase (black solid line) is known exactly 61 and serves as a validity check of our simulations (△). As for the FCSL ( and •), the exponential scaling law is found to be characterized by a parameter αFCSL = 0.473 ± 0.005 (dashed line). For the sake of completeness, the two cases where the initial pair of defects out of the FCSL were "2 in -2 out" (•) and "4 in / 4 out" ( ) were considered separately, but no pertinent difference was noticed. The y-axes are on a logarithmic scale, while the x-axes are on a a) linear and b) logarithmic scale to respectively emphasize the exponential scaling law and the long distance behavior. All data have been arbitrarily normalized to 1 at R = 1, r d being the shortest distance between charges. The error bars are smaller than the data symbols.

thus expected to come from the divergence-free contribution of the magnetization and, according to our results, to be smaller than in spin ice. Indeed, with χ being the susceptibility, the variance of the magnetization ∆M = (χT ) T →0 is 2.00(1) in spin ice [62][63][64] and 1.52(2) in the FCSL 22 with the same ratio as for the α's

∆M SI ∆M FCSL = 1.30(2) = α FCSL α SI (9)

IV. ENERGETIC COULOMB POTENTIAL FROM DIPOLAR INTERACTIONS

From now on, we shall include long-range dipolar interactions as given in Eq. (1).

Starting from a FCSL configuration and flipping one spin creates a pair of topological defects. But as opposed to spin ice, these defects can be categorized into two different kinds. On one hand, if a minority spin is flipped, then a pair of 4 in / 4 out topological defects is created, breaking the local divergence-free contribution of the Helmholtz decomposition (see Fig. 2). On the other hand, if a majority spin is flipped, then a pair of 2 in -2 out topological defects is created, breaking the local divergence-full contribution.

To study the nature of these topological defects, the averaged potential energy between them has been computed using the same worm algorithm as in section III. The total dipolar energy of the system is computed at each step of the defect diffusion and stored as a function of the distance between the two defects, until they eventually annihilate each other, forming a new FCSL configuration. The procedure is repeated for 10 4 worm updates, and averaged over 100 independent initial configurations, in order to obtain good statistics. Simulations were done for reasonably big system sizes (16 000 sites) in order to minimize finite size effects. Indeed, when the pair of defects is separated by half the system size, the influence of the mirror images due to periodic boundary conditions is not negligible. For a systematic analysis, the potential energy between two "2 in -2 out" defects (noted V hh ), between "4 in" and "4 out" defects (noted V mm ) and between "2 in -2 out" and "4 in / 4 out" defects (noted V hm ) were computed separately (see the three panels of Fig. 5).

Our results clearly confirm the emergent Coulomb potential between topological defects. It should be emphasized that this potential is directly coming from the dipolar energy between spins and should not be confused with the entropic potential considered in section III. Furthermore, as shown later in section VI, the FCSL degeneracy is weakly lifted by dipolar interactions. Even if small when compared to the total energy of the system, this degeneracy lift is extensive. It is thus remarkable that despite extensive fluctuations of the total dipolar energy of the system, especially for configurations with finite magnetization, the averaged interaction between a single pair of defects is quantitatively described by the Coulomb potential of Eq. ( 4).

Based on the dumbbell model 6 , it is possible to quantify these excitations analytically. When a pair of topological defects is created, the energy cost ∆E has three origins:

• the chemical potential of each defect, p h for "2 in -2 out" and p m for "4 in / 4 out";

• the Coulomb attraction between the newly cre- As computed in the Supplementary Information of Ref. [6] for spin ice, the dipolar interactions also contribute to the chemical potential of magnetic monopoles. When applied to the FCSL, the chemical potentials can be rewritten as

p h = p h + 8 3 1 + 2 3 D (10)
p m = p m -8 1 + 2 3 D, (11)
where p i={h,m} is the chemical potential of non-dipolar origin, such as from Hamiltonian (5) for example. As discussed in section II, the values of p i={h,m} are chosen such that the FCSL phase is stable. We now have all ingredients to estimate the energy cost for creating a pair of "2 in -2 out" defects (∆E hh ) or "4 in / 4 out" defects (∆E mm )

∆E hh = -2p h + V nn -2V M = -2p h + 16 3 -1 + M zb - 3 2 2 3 D = -2p h -4.73D, ( 12
) ∆E mm = -2p m + V nn + 2V M = -2p m + 16 + 8 3 2 3 (5 -2M zb ) D = -2p m + 19.75D. ( 13
)
In Monte Carlo simulations with dipolar interactions, corrections beyond the dumbbell model slightly modifies the numerical values of equations ( 12) and ( 13) and respectively give -4.34D and +19.70D; the comparison is excellent for ∆E mm and with a relative error of 8% for ∆E hh .

As for spin ice 6,65,66 , the emergent Coulomb potentials are thus rather natural in the context of the dumbbell model. But the motivation of this work is that the underlying magnetic order brings a new flavor to the problem. Once a pair of charges is created, additional spin flips separate the topological defects. Since the underlying long-range charge order takes place on a bipartite diamond lattice, with alternative positive and negative charges, a defect moving from one tetrahedron to the neighboring one alternatively becomes, for example, 2 in -2 out, 4 in, 2 in -2 out, 4 in ... carrying its topological and magnetic charge with it. This situation is analogous to kagome-ice physics 17,18 . Since "2 in -2 out" and "4 in" defects have a priori different chemical potentials (see Eqs. (10) and ( 11)), the hopping of topological defects from a tetrahedron to its neighbor alternatively costs and gains energy. Since the energy cost due to the chemical potential cancels out every two spin flips, this process does not confine the defects. Actually, this oscillating energy even disappears if we consider the simultaneous flipping of two spins successively pointing in and out, ensuring that the hopping of the quasi-particle remains on the same diamond sublattice. We shall refer to such dynamics as a double-spin motion, which is reminiscent of a dimer move on the diamond lattice and of the mobility of holons and spinons in frustrated Mott insulators on bipartite lattices [START_REF] Poilblanc | Mobile Holes in Frustrated Quantum Magnets and Itinerant Fermions on Frustrated Geometries[END_REF] . The creation of a pair of "2 in -2 out" and "4 in / 4 out" defects by a double-spin motion is depicted by green spins in the inset of Fig. 5.c), and costs an energy

∆E hm = -p h -p m + V nnn = -p h -p m + 4 + 16 √ 2 3 √ 3 D = -p h -p m + 8.35D (14)
where V nnn = -4D/3 is the next-nearest-neighbour Coulomb potential energy (see Eq. ( 4)). Please note that the Madelung term dissappears because the two defects with opposite charges sit on the same diamond sublattice. A Monte Carlo average with Ewald summation gives 8.67D, i.e. a difference of 4% with the numerical term of Eq. ( 14). When created alone, a pair of "4 in / 4 out" defects only carries ±1 single charges (see Fig. 5.b) as opposed to the ±2 double charges in spin ice. In the Helmholtz decomposition, this is because half of the zero-divergence has already been broken in the FCSL. However, when several topological defects are created, higher values of charges are possible. In the spin-ice model, the allowed magnetic charges are {-2, -1, 0, +1, +2}. Out of the FCSL, they are {-3, -2, -1, 0, +1, +2, +3}. If the longrange charge order of the FCSL imposes all up tetrahedra to be in a "3 out -1 in" spin configuration, then a "4 in" defect would carry a +3 magnetic charge on up tetrahedra. This means that, as opposed to spin ice, the magnetic charge carried by a topological defect is now defined in the context of its environment, and in particular by the charge-order symmetry breaking of the FCSL. The consequence is that even if "4 in / 4 out" defects always have the same sign -negative for 4 out, positive for 4 in -the value of their charge is not uniquely defined.

As for a 2 in -2 out defect, even the sign of its charge is unknown if considered alone. This is because a 2 in -2 out defect breaks the local divergence-full condition (see Fig. 2) and its sign is thus defined by the nature of the charge order that has been broken. If the charge order (i.e. divergence-full condition) covers the entire system apart from some dilute topological defects, then the sign of the defect is rather trivially given by the up/down na-FIG. 6. Schematic representation of the semiconductor analogy with the FCSL. The pseudo-magnetization ρ of Eq. ( 6) delimits the fraction of degrees of freedom respectively involved in the long-range ordered phase with divergence-full field (fully occupied valence band in yellow) and the Coulomb spin liquid with divergence-free field (empty conducting band). The divergence-free/full fields are explained in Fig. 2. A double-spin flip creates a gapped excitation (see Eq. ( 14) and inset of Fig. 5.c) taking the form of a pair of monopole and monopole hole which respectively break the local divergencefree and divergence-full fields. The gap is not visible in the figure because the monopole / monopole hole excitation conserves ρ. Further double-spin motions propagate these two quasi-particles within their respective conducting and valence bands, interacting only via an effective magnetic Coulomb potential (see Fig. 5.c). In 3 dimensions, this potential is non-confining and gives rise to fractionalization as in semiconductors. As for the "all in / all out" phase, the pseudomagnetization is saturated at ρ = 1. Topological defects thus necessarily take the form of monopole holes and are strongly confined since their diffusion requires the creation of further topological defects. In that sense, it is tempting to compare the pseudo-magnetization with the level of band filling: the hatched region represents an infinite gap because ρ cannot be higher than 1, and the "all in / all out" phase is analogue to an insulator of magnetic charges.

ture of the tetrahedron. However, if the system is covered by different domains of charge orderi.e. with different time-reversal symmetry breaking -then the immediate vicinity of a given 2 in -2 out defect is necessary to define the sign of its charge. The dynamics of the resulting domain walls, mediated by the co-existence with a Coulomb spin liquid and the emergence of topological magnetic charges would be a very interesting question to consider in the future.

V. SEMICONDUCTOR ANALOGY A. Monopole holes

As a summary, a pair of defects created out of the FCSL by a double-spin motion, or dimer move, is made of a "2 in -2 out" defect breaking the local divergence-full field, and a "4 in" or "4 out" defect breaking the local divergence-free field (see Fig. 2 and inset of Fig. 5.c). These defects carry a magnetic charge and interact via an emergent Coulomb potential (see Fig. 5.c). The fractionalization of topological defects in spin ice has already been noticed to carry some analogy with semiconductors 13,60 . In the FCSL phase, we can go one step further. In semiconductors, an electron hole is defined by the lack of an electron from what should have been a fully-occupied valence band. From the Helmholtz decomposition, a "2 in -2 out" defect can be understood as a monopole hole in the otherwise charge-ordered phase, while the "4 in / 4 out" defect is a monopole in the same sense as in spin ice, breaking the local divergence-free condition. Please note that each defect may carry either a positive or negative magnetic charge. In absence of defects, the long-range order (full-divergence) and spin liquid (zero-divergence) respectively play the role of the fully occupied valence band and empty conducting band of semiconductors, with gapped excitations to create a pair of "monopole / monopole hole" (see Eq. ( 14)). The 1/R potential ensures deconfinement of the topological defects and thus fractionalization of the excitations, even for monopole holes which break long-range order.

We should emphasize the importance of the doublespin motion in this analogy. Because of the time-reversal symmetry breaking of the long-range order, the divergence-full field of the FCSL is alternatively a source and a sink of fluxes on the bipartite diamond lattice. It means that if a negative magnetic charge is a monopole hole on up tetrahedra, it would be a monopole on down tetrahedra. Thus even if for a snapshot of the system at a given time t, the monopoles and monopole holes are properly defined, the nature of these defects is dynamically conserved only when hopping on the same diamond sublattice via double-spin motion. This is a direct consequence of the underlying dimer-model mapping. Of course here we are talking about Monte-Carlo dynamics. Real dynamics would be material dependent, but the system would likely take the form of a fluctuating ensemble of monopoles and monopole holes whose ratio would be determined by the temperature and the chemical potentials p i={h,m} , as if the doping was intrinsic rather than extrinsic.

An insightful way to appreciate the properties of the FCSL is to understand how it differs from the fully disordered ("2 in -2 out") and fully ordered ("all in / all out") phases. To do so, the pseudo-magnetization ρ is a useful observable, both globally for the entire system as defined in Eq. ( 6) and locally at the level of each tetrahedron (ρ tet = i∈tet σ i /4). In the "2 in -2 out" spin-ice ground state, the absence of magnetic order with ρ tet = 0 everywhere simply prevents the creation of monopole holes. This is why spin ice has often been discussed as a magnetic equivalent of an electrolyte 61,68 . Please note that both for spin ice 35,68,69 and the FCSL, only alternating currents of magnetic charges can subsist on long time scales. As for the "all in / all out" long-range order (ρ tet = 1 everywhere), this phase has recently raised much interest in the context of pyrochlore Iridates Nd 2 Ir 2 O 7 70 and Tb 2 Ir 2 O 7 71 . In comparison with the FCSL, "all in / all out" order means that the conducting band has been fully filled by divergence-full monopoles, making the pseudo-magnetization an effective measure of the level of band filling, as illustrated in Fig. 6. Interestingly, defects out of the "all in / all out" phase are necessarily confined monopole holes with no empty band available since ρ cannot be higher than 1. In that sense, the "all in / all out" phase is similar to an insulator of magnetic charges with an infinite gap. A possible follow-up of this analogy would be to study how this "insulator" could couple with a nearby Coulomb phase, such as for example in FeF 3 72 or possibly with the presence of two interpenetrating magnetic pyrochlore lattices such as in pyrochlore Iridates.

The framework we have developed here provides a working example of the kind of emergent phenomena arising from the co-existence of magnetic order and spin liquids. The semiconductor analogy is able to catch a surprising number of characteristic features of the FCSL. There are of course limitations to this analogy, which we believe are good places to look for exotic physics, such as the possibility of domain walls separating domains with opposite time-reversal symmetry or the inherent topological nature of the FCSL which would probably lead to phase transitions beyond the standard Landau-Ginzburg-Wilson paradigm [73][74][75][76][77] .

B. Beyond pyrochlores

So far, we have focused on the pyrochlore lattice in order to develop a comprehensive background for our theory. But the present analogy between semiconductor physics and a partially ordered Coulomb phase is not restricted to the pyrochlore symmetry. The microscopic mechanism able to stabilize a FCSL would always remain a model-dependent issue, but conceptually, our theory can be easily extended to most vertex models supporting a Coulomb phase. A straightforward ramification would be the bipartite hexagonal (Ih) ice lattice, which is the standard form of water ice at ambient pressure 48 and of recent interest for quantum proton dynamics [78][79][80][81] . At lower dimensions, the same analogy can be made on the checkerboard lattice but with the caveats that entropic correlations are logarithmic in two dimensions, FIG. 7. Vertices on the cubic lattice have 6-fold connectivity as opposed to the tetragonal geometry of the diamond lattice. This higher connectivity gives rise to two kinds of FCSL, i.e. with co-existence of long-range order (ρ > 0) and Coulomb phase: the hard-core dimer model at ρ = 2/3 (dimers are in green) and the fully packed loop model at ρ = 1/3 (the loop is a thick red line). Upon monopoles doping (i.e. increasing ρ), it is thus possible to visit two phases similar to semiconductors before encountering the "insulating" all in / all out long-range order at ρ = 1.

while magnetic dipolar interactions are known to lift the "2 in -2 out" degeneracy at the vertex level, a common problem in artificial spin ice 82,83 . The latter issue is absent in kagome ice, and so is the equivalent "2 in -2 out" phase since a proper definition of the divergence-free field on kagome requires time-reversal symmetry breaking. Interestingly, a FCSL is known to be stable over a finite temperature window in presence of dipolar interactions on kagome 17,18 , or equivalently on the [111] magnetization plateau of spin ice materials 7,84 .

As illustrated by kagome, the connectivity of the vertex does not need to be tetragonal. Let us thus consider the bipartite cubic lattice which has been recently studied in the context of unconventional phase transitions in dimer/monomer and loop models 85,86 and emergent loop models in perovskite oxynitrides 87 .We do not con-sider long-range dipolar interactions here, and thus the defects do not carry any effective magnetic charges. We impose that each bond of the cubic lattice carries a magnetic flux represented by an arrow towards one of the two neighbouring vertices. As opposed to the diamond lattice, each vertex has 6 legs (or arrows). By imposing all vertices of a given sublattice to have the same number of entering fluxes (inwards arrows), one can define various Coulomb phases as depicted in Fig. 7:

• 3 entering fluxes: there is no long-range order (ρ = 0).

• 2 entering fluxes: there is partial order (ρ = 1/3) and the phase is mapped onto the fully packed loop model without loop crossing. Fluctuations of the loops represent fluctuations of the co-existing Coulomb spin liquid.

• 1 entering flux: there is partial order (ρ = 2/3) and the phase is mapped onto the hard-core dimer model.

• 0 entering flux: the system is fully ordered (ρ = 1).

On the pyrochlore lattice, when starting from the FCSL (ρ = 1/2) and filling the conducting band with "standard" monopoles, one would immediately reach the insulating "all in / all out" order (ρ = 1). On the cubic lattice on the other hand, starting from ρ = 1/3 and filling the conducting band, one would first reach another FCSL with ρ = 2/3, offering the tempting analogy with multi-band physics.

VI. GROUND STATE OF THE DIPOLAR FCSL: MAXIMIZATION OF FLIPPABLE PLAQUETTES

Even if the potential energy between a pair of defects obeys the Coulomb law on average (see Fig. 5), the ensemble of configurations of the FCSL remains quasidegenerate. It means that at low enough temperature, the system should ultimately order.

As shown in Fig. 8, Monte Carlo simulations within the FCSL ensemble (without any topological defect) confirm the nature of the ground state found in Ref. [21], which we identify as the so-called R-states depicted in Fig. 9. The phase transition is of first order and occurs at T c /D = 0.166±0.003. Coming from high temperature, the build up of correlations due to the R-states is clearly visible in the Non-Spin-Flip channel of the structure factor (see Fig. 10), while the inherent divergence-full longrange order of the FCSL is signaled by well-defined Bragg peaks (black dots) in the Spin-Flip channel 22 . Remarkably, the pinch points due to the divergence-free Coulomb phase persist even just above the transition temperature.

To find the R-state here is quite interesting. Indeed, the R-state has been first studied to explain the magnetisation plateau of HgCr The hysteresis of the transition makes it difficult to precisely estimate the transition temperature. This is why in order to show convergence of the simulations, two kinds of equilibration processes were used. The system was either quenched into one of the the R-states (•) or slowly annealed from high temperature (△), which respectively over-and underestimate the transition temperature. The MC time tMC of the x-axis represents the number of worm updates which are accepted/rejected based on a Metropolis argument accounting for the dipolar interactions. For both equilibration processes, the simulations were equilibrated at the temperature of measurements during tMC /10. The results are given for different system sizes, showing convergence of the simulations towards Tc/D = 0.166 ± 0.003.

been shown to be the eight-fold degenerate ground state of the quantum dimer model on the diamond lattice in the limit of large negative potential energy [24][25][26] .

In the quantum dimer language, a negative potential energy favor configurations with flippable plaquettes, i.e. where rings of 3 dimers can dynamically resonate around a hexagon of 6 sites on the diamond lattice. As explained in Fig. 1, the FCSL is the classical counterpart of this quantum dimer model (modulo time-reversal symmetry), where the minority spins are the dimers. In our spin language, a flippable plaquette is represented by thick lines in Fig. 9; flipping the 6 spins around this hexagon is the smallest fluctuation of the Coulomb spin liquid keeping the long-range charge order intact, and is an excitation out of the R-states.

Our physical understanding of this analogy is that dipolar interactions, by definition, enforce the closing of magnetic flux lines. As a textbook example, dipolar interactions are the reason for Weiss domains in standard ferromagnets. Since these "flippable plaquettes" are the smallest ensemble of spins able to close a magnetic flux line, it is understandable that dipolar interactions may want to maximize the number of such plaquettes. At the classical level, this is similar to a negative potential energy in the quantum dimer model, which is why the R-states are found as the ground states of both models.

But this result raises a new question. The ground state of spin ice also has a quantum counterpart, which has been studied in terms of potential and kinetic energies of flippable plaquettes [1][2][3][4] . But while dipolar spin ice orders in a q = (001) antiferromagnetic ground state 38 , its quantum counterpart orders into the so-called "squiggle" state 88 for large negative potential energy. We propose that the reason of this difference comes from the longrange nature of the dipolar interactions, which favor the closing of the flux lines for all length scales. In particular, a finite macroscopic magnetization is strongly hindered by dipolar interactions. This is not the case for the potential term of the quantum model which is local. Hence the finite magnetization of the squiggle state along the (001) axis 88 is forbidden by dipolar interactions, but not in the quantum model. Since the R-states are antiferromagnetic, we do not encounter this problem in the ordering of the FCSL.

VII. CONCLUSION

The co-existence of partial magnetic order and spin fluctuations has been a recurrent experimental and often controversial observation over the past years, especially in rare-earth pyrochlores Yb 2 Ti 2 O 7 89-92 , Tb 2 Ti 2 O 7 93,94 , Er 2 Sn 2 O 7 [95][96][97][98] and Vesignieite kagome [99][100][101] . Here we have developed a theoretical framework to characterize the co-existence of long-range order and a Coulomb gauge field which we have named a Fragmented Coulomb Spin Liquid.

As for excitations out of the spin-ice ground state, topological defects of the FCSL interact via effective entropic and energetic Coulomb interactions, the latter one being due magnetic dipolar interactions between spins. Following the Helmholtz decomposition, topological defects can be categorized into two kinds: monopoles which break the local divergence-free field as in spin ice, and monopole holes which break the divergence-full field and are thus directly related to the co-existence with longrange order. Our work i) shows that partial long-range order does not necessarily prevent deconfinement of fractionalized excitations, which would be consistent with persistent dynamics, and ii) shed a new light on the monopole picture in analogy with semiconductors, as a direct consequence of the underlying long-range order. Being based on the simple decomposition of degrees of freedom between an ordered contribution and a fluctuating one, we expect this framework to be relevant to a wide range of systems (artificial spin ice, magnetisation plateaux...).

Finally, we have identified that dipolar interactions on the classical FCSL support the same ground state as large negative potential energy on the equivalent quantum dimer model 21,[24][25][26] . This offers a natural classical mechanism to maximize the number of flippable plaquettes, which can be expected to be a general tendency of dipolar interactions, but not a systematic outcome; among other points, it depends on the longdistance geometric constraints. Future directions of work would be to consider Heisenberg spins instead of Ising, which would for example allow for gapless excitations. At the classical level, the decomposition of degrees of freedom would be straightforward in terms of irreducible representations 98 : the divergencefull and divergence-free contributions would respectively correspond to A 2 and T 1 symmetries. Also, instead of one long-range antiferromagnetic phase, we could consider the influence of frozen disorder 102-104 , or nematic ) and eSF vectors respectively. eSF is defined by forming an orthogonal basis with eNSF and the wave-vector q. Both temperatures are above the transition temperature of Tc/D = 0.166. Black dots in the left panels are Bragg peaks signaling the long-range charge order of the FCSL 22 . Despite the build-up of R-states correlations in the NSF channel, pinch points are still visible in the total scattering. phases where quadrupolar decomposition might be necessary.

FIG. 5 .• 8 3 2 3 D

 5833 FIG.5. Averaged potential energy between a) two "2 in -2 out" defects (•), b) "4 in" and "4 out" defects (△), c) "2 in -2 out" and "4 in / 4 out" defects ( ), as illustrated in the inset of each panel, obtained from Monte Carlo simulations with dipolar interactions. "2 in -2 out" tetrahedra are colorless, while the "4 in / 4 out" tetrahedra are represented by thick lines. The green spins in panel c) represent the double-spin motion responsible for the creation of the "2 in -2 out" and "4 in / 4 out" defects. To compare the simulations with Eq. (4) (dashed lines), the data have been shifted by a reference energy, serving as fitting parameter, whose origin is discussed in Eqs.(12-14). There are no multiplicative prefactors, meaning that all defects carry a unit of magnetic charge. The dipolar energy scale D is set to 1. The error bars, computed from the results of 100 independent samples, are smaller than the data symbols.

FIG. 8 .

 8 FIG.8. At very low energy scale, dipolar interactions ultimately lift the FCSL degeneracy at Tc/D = 0.166 ± 0.003. a) The first-order nature of the transition is clearly visible in the specific heat C h (•) and in the discontinuity of the R-state order parameter ( ). The data are obtained from Monte Carlo (MC) simulations within the FCSL ensemble, averaged over 10 7 worm updates and 10 independent initial configurations, for a system size of 3456 pyrochlore sites. b) The hysteresis of the transition makes it difficult to precisely estimate the transition temperature. This is why in order to show convergence of the simulations, two kinds of equilibration processes were used. The system was either quenched into one of the the R-states (•) or slowly annealed from high temperature (△), which respectively over-and underestimate the transition temperature. The MC time tMC of the x-axis represents the number of worm updates which are accepted/rejected based on a Metropolis argument accounting for the dipolar interactions. For both equilibration processes, the simulations were equilibrated at the temperature of measurements during tMC /10. The results are given for different system sizes, showing convergence of the simulations towards Tc/D = 0.166 ± 0.003.

FIG. 9 .

 9 FIG. 9. Unit cell of the eight-fold degenerate ground state of the FCSL with dipolar interactions, as observed in Ref.[21], which we identify as the R-states. Within a unit cell, the magnetization of the four down tetrahedra point in all four different local (111) directions. The R-state is known for maximizing the number of flippable plaquettes within the FCSL. The minority spins are colored in green and the flippable plaquette is shown by thick lines.

FIG. 10 .

 10 FIG. 10. Structure factor of the FCSL with dipolar interactions obtained from Monte Carlo simulations at T /D = 2 (top) and T /D = 0.225 (bottom), showing the Spin-Flip (left), Non-Spin-Flip (middle) and total (right) scattering. The NSF and SF channels measure the spin components along the eNSF = (1 10) and eSF vectors respectively. eSF is defined by forming an orthogonal basis with eNSF and the wave-vector q. Both temperatures are above the transition temperature of Tc/D = 0.166. Black dots in the left panels are Bragg peaks signaling the long-range charge order of the FCSL22 . Despite the build-up of R-states correlations in the NSF channel, pinch points are still visible in the total scattering.

ACKNOWLEDGMENTS

The author is thankful to Marion Brooks-Bartlett, Simon Banks, Laurent de Forges de Parny and Peter Holdsworth for collaborations on related topics and to Owen Benton, Andrew Smerald and Roderich Moessner for insightful comments on the manuscript. This work was supported by the Okinawa Institute of Science and Technology Graduate University.