# Competing Spin Liquids and Hidden Spin-Nematic Order in Spin Ice with Frustrated Transverse Exchange

Mathieu Taillefumier, Owen Benton, Han Yan, Ludovic D.C. Jaubert, Nic Shannon

► To cite this version:

Mathieu Taillefumier, Owen Benton, Han Yan, Ludovic D.C. Jaubert, Nic Shannon. Competing Spin Liquids and Hidden Spin-Nematic Order in Spin Ice with Frustrated Transverse Exchange. Physical Review X, 2017, 7, pp.041057. 10.1103/PhysRevX.7.041057. hal-01549587

HAL Id: hal-01549587

https://hal.science/hal-01549587

Submitted on 28 Jun 2017

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Frustrating quantum spin ice : a tale of three spin liquids, and hidden order

Mathieu Taillefumier,¹ Owen Benton,² Han Yan,¹ Ludovic Jaubert,³ and Nic Shannon¹

¹Okinawa Institute of Science and Technology Graduate University, Onna-son, Okinawa 904-0495, Japan

²RIKEN Center for Emergent Matter Science (CEMS), Wako, Saitama, 351-0198, Japan

³CNRS, Université de Bordeaux, LOMA, UMR 5798, 33400 Talence, France

Frustration in magnetic interactions often gives rise to disordered ground states with subtle and beautiful properties. The spin ices Ho₂Ti₂O₇ and Dy₂Ti₂O₇ exemplify this phenomenon, displaying a classical spin liquid state, with fractionalized magnetic-monopole excitations. Recently there has been great interest in closely-related “quantum spin ice” materials, following the realization that anisotropic exchange interactions could convert spin ice into a massively-entangled, quantum, spin liquid, where magnetic monopoles become the charges of an emergent quantum electrodynamics. Here we show that by allowing for frustration in the transverse interactions, on top of the inherent frustration of spin ice, an even richer scenario can be realized. Using large-scale numerical simulations of the simplest model of a quantum spin ice, the XXZ model on the pyrochlore lattice, we find that frustrated transverse exchange opens the door to not one, but three distinct types of spin liquid, as well as a phase with hidden, spin-nematic, order. We explore the experimental signatures of each of these new phases, making explicit predictions for inelastic neutron scattering. These results show an intriguing similarity to experiments on a range of pyrochlore oxides.

## I. INTRODUCTION

The search for spin liquids — disordered phases of magnets which support entirely new forms of magnetic excitation — has become one of the defining themes of modern condensed-matter physics [1, 2]. In this context, the pyrochlore lattice, a corner-sharing network of tetrahedra found in a wide range of naturally-occurring minerals, has proved an amazing gift to science. Pyrochlore magnets play host to a variety of unconventional forms of magnetic order, and include systems which have not been observed to order at any temperature [3]. Perhaps the most celebrated of these is the “spin ice” found in the Ising magnets Ho₂Ti₂O₇ and Dy₂Ti₂O₇ [4]; a classical spin liquid, described by an emergent U(1) lattice gauge theory with magnetic monopole excitations [5].

As the understanding of spin ice has grown, so more attention has been given to the role of quantum effects. These are of particular relevance where a spin-ice arises through anisotropic exchange interactions in a pyrochlore magnet [6–11], and have the potential to convert classical spin ice into a massively-entangled, quantum, spin liquid, described by an emergent U(1) quantum electrodynamics [11–21]. At the same time, great progress has been made in synthesizing and characterizing magnetic pyrochlore oxides. As well as revealing a number of candidates for quantum spin-ice behaviour [22–29], these experiments have turned up many unusual and unexpected magnetic states in systems with strongly anisotropic exchange [30–39].

The main message of this Article is that even the simplest model of a quantum spin ice — the XXZ model on a pyrochlore lattice — has far more to offer than spin ice alone. Concentrating on frustrated interactions, in the classical limit, accessible to large-scale simulation, we find not one, but three, distinct types of spin liquid, each of which is described by a different lattice gauge theory [Fig. 1]. As an added bonus, one of these spin

liquids undergoes a phase transition at low temperatures into a state with hidden, spin-nematic, order. Explicit predictions are made for the experimental signatures of each of these different spin liquids [Fig. 2], and molecular dynamics simulations are used to characterise the gapless collective excitations of the spin-nematic phase [Fig. 5]. We find that the portrait which emerges has striking similarities with the behavior of a number of pyrochlore materials.

The simplest model able to capture quantum effects in a spin ice is the XXZ model on the pyrochlore lattice

$$ \mathcal{H}_{\text{QSI}} = \sum_{\langle ij \rangle} J_{zz} S_i^z S_j^z - J_{\pm} (S_i^+ S_j^- + S_i^- S_j^+) \quad (1) $$

where $S_i = (S_i^x, S_i^y, S_i^z)$ is a (pseudo)spin-half operator describing the two states of the lowest energy doublet of a magnetic ion subject to a strong crystal electric field (CEF). The symmetry of the lattice requires that the quantization axis of each spin (here, $S_i^z$) lies on a local [111] axis, as defined in Appendix A.

Ising interactions, $J_{zz} > 0$, favor states obeying the “ice rules” in which two spins point into, and two spins point out of, each tetrahedron on the lattice. The transverse term, $J_{\pm}$, introduces dynamics about these spin-ice configurations and, for larger values of $J_{\pm}/J_{zz}$, can drive the system into a state with easy-plane order [16, 17, 40–43]. The physical meaning of this easy-plane order depends on the nature of the magnetic ion. For Kramers ions like Yb³⁺ and Er³⁺ all components of S relate to a magnetic dipole moment [10], and the ordered phase is an easy-plane antiferromagnet. However for non-Kramers ions such as Pr³⁺ and Tb³⁺ [8, 38], or “dipolar-octupolar” Kramers ions like Nd³⁺ or Ce³⁺ [44], the easy-plane order may have quadrupolar (octupolar) character. In what follows, we consider explicitly the case of Kramers ions.

For $J_{\pm} > 0$, $\mathcal{H}_{\text{QSI}}$ [Eq. (1)] is unfrustrated, in the sense that it is free of sign problems in Quantum Monte Carlo
<figure><img src="image_2.png" /><figcaption>Figure 1. (a) Finite-temperature phase diagram of the XXZ model on a pyrochlore lattice, $\mathcal{H}_{\text{QSL}}$ [Eq. (1)]. The model possesses three distinct spin-liquid phases (SL), as well as ground states with easy-plane antiferromagnetic ($AF_{\perp}$) and spin-nematic ($Q_{\perp}$) order. Associated crossover temperatures (phase transitions) are indicated with dashed (solid) lines. For $J_{\pm}/J_{zz} = -1/2$ (dash-dotted line), the model is thermodynamically equivalent to the Heisenberg antiferromagnet on a pyrochlore lattice. Results are taken from classical Monte Carlo simulation of a cubic cluster of $N = 8192$ spins. (b) Representative configuration of quadrupoles in the $Q_{\perp}$ phase with nematic order (c) “Two in, two out” configuration of spins in the spin ice regime (SLI). (d) Representative configuration of spin dipoles in the ordered $AF_{\perp}$ phase.</figcaption></figure>

(QMC) simulation. In this case, the phase diagram is already well-established [13, 20, 21]. For $J_{\pm}/J_{zz} \lesssim 0.05$, QMC simulations find a crossover from a conventional paramagnet into a classical spin-liquid (spin ice) at a temperature $T^*/J_{zz} \sim 0.2$, and a second crossover into a quantum spin liquid (QSL) at a much lower temperature $T^*_{QSL}/J_{zz} \sim (J_{\pm}/J_{zz})^3$. In the low temperature quantum spin liquid regime, the magnetic monopoles of classical spin ice become dynamic, fractional, spin excitations (spinons), while the spectrum of the model also includes gapless photons [12, 15]. For $J_{\pm}/J_{zz} \gtrsim 0.05$, the U(1) QSL gives way to easy-plane antiferromagnetic order ($AF_{\perp}$), in which spins lie in the plane perpendicular to the local $S^z$-axis [13, 20, 21].

Very little is known about the properties of $\mathcal{H}_{\text{QSL}}$ for $J_{\pm} < 0$ [9, 17, 38]. On perturbative grounds, it is expected that the ground state for $|J_{\pm}|/J_{zz} \ll 1$ will also be a U(1) QSL [12], albeit one with a modified spinon dis-

persion [17, 45]. Gauge Mean-field calculations suggest that this QSL persists over a broad range of parameters, $-4.13 \lesssim J_{\pm}/J_{zz} < 0$ [17]. But the nature of competing ordered — or disordered — phases for $J_{\pm} < 0$ remains an open question.

There are many reasons to believe that the properties of the quantum spin ice model, $\mathcal{H}_{\text{QSL}}$ [Eq. (1)] for frustrated coupling $J_{\pm} < 0$, could be even richer than for $J_{\pm} > 0$. In particular, for $J_{\pm}/J_{zz} = -\frac{1}{2}$, $\mathcal{H}_{\text{QSL}}$ [Eq. (1)] is equivalent (up to a site-dependent spin-rotation), to the Heisenberg antiferromagnet (HAF) on a pyrochlore lattice. Like spin ice, the HAF is known to support a classical spin liquid [46–49], and it has also been argued to support a QSL ground state [50–53]. And, crucially, both the classical and quantum spin liquids in the HAF have a qualitatively different character from those found in spin ice. This sets up a competition between two different kinds of spin liquid, namely spin ice for $J_{\pm}/J_{zz} \approx 0$,
and a state homologous to the HAF for $J_{\pm}/J_{zz} \approx -\frac{1}{2}$. It also opens the door for yet more novel magnetic phases for $J_{\pm}/J_{zz} < -\frac{1}{2}$.

## II. PHASE DIAGRAM DETERMINED BY CLASSICAL MONTE CARLO SIMULATIONS

Since the quantum spin ice model, $\mathcal{H}_{QSI}$ [Eq. (1)], is inaccessible to QMC for $J_{\pm} < 0$, we instead study its finite-temperature properties using classical Monte Carlo (MC) simulation — the results are summarised in the phase diagram Fig. 1. For $J_{\pm} > 0$, this phase diagram is very similar to that previously found in QMC simulations [13, 20, 21] — at a qualitative level, the only significant difference is the absence of a QSL below $T_{QSL}^*/J_{zz} \sim (J_{\pm}/J_{zz})^3 \lesssim 0.005$. At a quantitative level, we find changes in numerical values of the crossover temperature associated with the spin ice regime, $T_1^*$, and the parameters for zero-temperature boundary between SL I and AF$_{\perp}$. These changes can be ascribed to the fact that the magnetic monopoles (spinons) are not quantized in classical simulations and do not develop phase coherence [55]. Further details of classical MC simulations for $J_{\pm} > 0$ will be presented elsewhere [56].

We now turn to the frustrated case, $J_{\pm} < 0$. At low temperatures, spin-ice correlations persist up to $J_{\pm}/J_{zz} = -\frac{1}{2}$ [9, 38], as illustrated in Fig. 2a. Upon reaching $J_{\pm}/J_{zz} = -\frac{1}{2}$ the system becomes thermodynamically equivalent to a HAF. This high-symmetry point gives rise to a new form of spin liquid at finite temperature, labelled SL III in Fig. 1. Once again, this spin liquid has algebraic correlations, as shown in Fig. 2c, but with qualitatively different character from spin ice [Fig. 2a]. These correlations persist up to a crossover temperature $T_3^*$ associated with the Curie-law crossover (CLC) in the magnetic susceptibility [57].

While the correlations measured in the equal-time structure factor $S(q)$ are also different from those found in the HAF [47, 48, 58], the two models are equivalent up to a local coordinate transformation. And, by analogy with earlier work on the HAF [48, 49, 59], the spin liquid SL III can be described by a $U(1) \times U(1) \times U(1)$ gauge theory. We refer to this state below as the pseudo-Heisenberg antiferromagnet (pHAF).

The situation for $J_{\pm}/J_{zz} < -\frac{1}{2}$ is even more interesting. Below a second crossover scale, $T_2^* < T_3^*$, identifiable by a reduction in the fluctuations of the z-components of the spins [see Appendix B], the spin liquid SL III (pHAF), gives way to a new spin liquid, labelled SL II in Fig. 1. Spin correlations in this phase have algebraic character, with pinch-points in $S(q)$ [Fig. 2b]. However these correlations are qualitatively different from those in either spin ice [Fig. 2a], or the pHAF [Fig. 2c]. At a still lower temperature, $T_Q < T_2^*$, the system undergoes thermodynamic phase transition, marked by a clear anomaly in the specific heat. None the less, this phase transition does not give rise to any magnetic Bragg peaks in $S(q)$ and,

at least as far as dipolar spin correlations are concerned, the system remains disordered.

While the new phase for $T < T_Q$ — labelled $Q_{\perp}$ in Fig. 1 — does not exhibit any conventional magnetic order, it does have a hidden, spin-nematic, order. The associated wavevector is $\mathbf{q} = 0$, and the relevant order parameter is,

$$ \mathbf{Q}_{\perp} = \left( \frac{\mathbf{S}^{x2} - \mathbf{S}^{y2}}{2\mathbf{S}^{x}\mathbf{S}^{y}} \right), \qquad (2) $$

where $S_i = (S_i^x, S_i^y, S_i^z)$ is written in the local frame of site $i$, defined in Appendix A. This order parameter quantifies the spontaneous breaking of the global $U(1)$ symmetry by the singling out of one axis within the local $xy$-plane. $Q_{\perp}$ transforms with the E–irrep of the tetrahedral point-group $T_d$, and is formally identical to the order parameter for the spin–nematic phases found in some frustrated magnets in applied magnetic field [60–62]. In common with these systems, the associated Landau theory

$$ F_{Q_{\perp}} = a_2(T) \mathbf{Q}_{\perp}^2 + a_4 \mathbf{Q}_{\perp}^4 + \dots, \qquad (3) $$

lacks a cubic term, and therefore permits a continuous phase transition. Simulations suggest that the phase transition at $T = T_Q$ is indeed continuous for $J_{\pm}/J_{zz} \lesssim -\frac{1}{2}$, becoming first-order approaching the high-symmetry point $J_{\pm}/J_{zz} \to -\frac{1}{2}$. Further details of the thermodynamics of this transition are given in Appendix B.

## III. THEORY OF THE SPIN LIQUID SL II

Spin correlations in spin ice (SL I) can be described using a $U(1)$ lattice gauge theory [5, 49, 63], which gives rise to characteristic “pinch-points” in the spin structure factor $S(q)$ [Fig. 2a]. Meanwhile spin correlations in the HAF — and by extension SL III (pHAF) — can be described using a $U(1) \times U(1) \times U(1)$ gauge theory [46–49]. This has qualitatively different pinch-points, as illustrated in Fig. 2c. It is clear that the correlations of SL II [Fig. 2b] are very different from either SL I or SL III. None the less, the presence of pinch points suggests that SL II, too, may be described by some form of gauge theory.

We can develop a field-theory for the spin liquid SL II by applying the methods developed in [43, 64] to treat a different kind of spin liquid. The starting point of this approach is to recast the spins $S_i$ in $\mathcal{H}_{QSI}$ [Eq. (1)] in terms of five order-parameter fields

$$ \{\mathbf{m}_{\lambda}\} = \{m_{A_2}, \mathbf{m}_{E}, \mathbf{m}_{T_1\text{ice}}, \mathbf{m}_{T_1\text{planar}}, \mathbf{m}_{T_2}\} \qquad (4) $$

defined on each tetrahedron $\mathbf{r}$. These objects $\mathbf{m}_{\lambda}(\mathbf{r})$ describe the different kinds of four-sublattice magnetic order consistent with the point group symmetry of the pyrochlore lattice. Definitions of each field $\mathbf{m}_{\lambda}$ in terms of the spins $S_i$ are given in Appendix D.

The most general exchange Hamiltonian on the pyrochlore lattice can be transcribed exactly in terms of $\mathbf{m}_{\lambda}$
<figure><img src="image_3.png" /><figcaption>Figure 2. Comparison of correlations in spin–liquid (SL) and spin–nematic ($Q_{\perp}$) phases. a) Spin ice (SL I), showing “pinch points” indicative of algebraic spin correlations. Left half-panel: equal-time structure factor $S(q)$, as measured in unpolarised neutron-scattering experiments. Right half-panel: $S(q)$ resolved into spin-flip (SF, top) and non spin-flip (NSF, bottom) components, as measured in polarised neutron-scattering experiments (cf. [54]). Definitions of each structure factor are given in Appendix C. b) Phase with spin–nematic order ($Q_{\perp}$), showing the absence of magnetic Bragg peaks, and persistence of algebraic spin correlations, which is also characteristic of the spin liquid SL II. c) Spin liquid SL III, for parameters thermodynamically equivalent to a Heisenberg antiferromagnet, showing algebraic spin correlations distinct from those in spin ice (SL I) or SL II. Results are taken from Monte Carlo simulations of $\mathcal{H}_{QSI}$ [Eq. (1)], for a cubic cluster of $N = 8192$ spins.</figcaption></figure>

[42]. This greatly simplifies the determination of classical ground states and, where classical ground states form an extensive manifold, one can use this approach to determine the local constraints which control the resulting spin–liquid [43, 64]. In the case of SL II, for $T \to 0$, we have

$$m_{A_2}(\mathbf{r}) = 0, \quad m_{E}(\mathbf{r}) = 0, \quad m_{T_1\text{ice}}(\mathbf{r}) = 0. \quad \forall \mathbf{r} \quad (5)$$

The spin fluctuations at low temperature are thus dominated by the fluctuations of the remaining fields $m_{T_2}(\mathbf{r})$ and $m_{T_1\text{planar}}(\mathbf{r})$. These fields have significance as the order-parameters of the competing four-sublattice magnetic orders which would be induced by the symmetry–allowed perturbation

$$\delta\mathcal{H}_{\pm\pm} = \sum_{\langle ij \rangle} J_{\pm\pm} [\gamma_{ij}S_i^+S_j^+ + \gamma_{ij}^*S_i^-S_j^-], \quad (6)$$

where $\gamma_{ij}$ are complex phase factors arising from the change in coordinate frame between different lattice sites [6, 8–10, 40, 65]. For this reason, the spin–liquid SL II falls very naturally into the “multiple–phase competition” scenario for pyrochlore magnets [42, 43, 66, 67].

In Fig. 3, we show the classical ground–state phase diagram of anisotropic exchange model

$$\mathcal{H}_{\text{ex}} = \mathcal{H}_{\text{QSI}} + \delta\mathcal{H}_{\pm\pm}. \quad (7)$$

This contains three distinct regions of 4–sublattice order: the easy–plane ordered phases described by the fields $m_E$ (denoted $AF_{\perp}$ in Fig. 1), $m_{T_1\text{planar}}$, and $m_{T_2}$ (Palmer–Chalker state [68]). These border a region of spin ice (denoted SL I in Fig. 1), dominated by fluctuations of $m_{T_1\text{ice}}$. We note that a closely–related phase diagram has been derived for non–Kramers ions [9, 38]; in this

case easy–plane order must be interpreted in terms of the quadrupole moment of the magnetic ion.

The non-trivial correlations in the spin–liquid SL II arise from the fact that neighbouring tetrahedra share a spin, so that the fields $m_{\lambda}(\mathbf{r})$ on neighbouring tetrahedra are not independent of one another. This point, combined with Eq. (5), imposes spatial constraints on the fluctuations of $m_{T_2}(\mathbf{r})$ and $m_{T_1\text{planar}}(\mathbf{r})$. After coarse graining to extract the long wavelength physics these constraints may be written in terms of two, independent, vector fluxes

$$\begin{align}
\mathbf{B}_1 &= \frac{1}{2} (2m_{T_1\text{planar}}^x, -\sqrt{3}m_{T_2}^y - m_{T_1\text{planar}}^y, \nonumber \\
&\phantom{=} \sqrt{3}m_{T_2}^z - m_{T_1\text{planar}}^z) \nonumber \\
\mathbf{B}_2 &= \frac{1}{2} (2m_{T_1\text{planar}}^x, -m_{T_2}^y + \sqrt{3}m_{T_1\text{planar}}^y, \nonumber \\
&\phantom{=} -m_{T_2}^z - \sqrt{3}m_{T_1\text{planar}}^z), \tag{8}
\end{align}$$

which each separately obey their own Gauss’ law

$$\nabla \cdot \mathbf{B}_1 = 0, \quad \nabla \cdot \mathbf{B}_2 = 0. \quad (9)$$

We can therefore write

$$\mathbf{B}_1 = \nabla \times \mathbf{A}_1, \quad \mathbf{B}_2 = \nabla \times \mathbf{A}_2 \quad (10)$$

and the theory has two, independent, U(1) gauge degrees of freedom.

The free energy associated with the fluctuations of these fields is of entropic origin [59]. The only choice of Gaussian free-energy consistent with both the point group symmetry and the U(1) symmetry of $\mathcal{H}_{QSI}$ is

$$\begin{align}
\mathcal{F}_{\text{SL II}} &= \frac{T}{V} \int d^3\mathbf{r} \, \lambda (\mathbf{B}_1^2 + \mathbf{B}_2^2) \nonumber \\
&= \frac{T}{V} \int d^3\mathbf{r} \, \lambda [(\nabla \times \mathbf{A}_1)^2 + (\nabla \times \mathbf{A}_2)^2] \tag{11}
\end{align}$$
<figure><img src="image_4.png" /><figcaption>Figure 3. (Color online). Classical ground-state phase diagram of the anisotropic exchange model $\mathcal{H}_{ex}$ [Eq. (7)], for $J_{zz} > 0$. Different phases are labelled in terms of the irreps of the tetrahedral symmetry group, $T_d$ [Eq. (4)], as described in [42]. The minimal model of a quantum spin ice $\mathcal{H}_{QSI}$ [Eq. (1)] exists on the line $J_{\pm\pm} = 0$ — for $J_{\pm} < -\frac{1}{2}$ (white line), two phases with 4-sublattice easy-plane order meet, and the resulting enlarged ground-state manifold gives rise to the spin liquid SL II, and spin-nematic phase $Q_{\perp}$. A closely-related mean-field phase diagram for non-Kramers ions is given in [9, 38].</figcaption></figure>

where the coefficient $\lambda$ can be determined through fits to simulation, or a large-N expansion [48, 64].

It follows from the existence of the conserved fluxes $\mathbf{B}_1$ and $\mathbf{B}_2$ and the free-energy Eq. (11) that SL II is a Coulomb phase with algebraic correlations [59]. The validity of this description is demonstrated in Fig. 4 where we compare analytic calculations of the flux structure factor

$$S_{\mathbf{B}_{\mu}}^{\alpha\beta}(\mathbf{q}) = \langle B_{\mu}^{\alpha}(-\mathbf{q}) B_{\mu}^{\beta}(\mathbf{q}) \rangle \quad (12)$$

based on Eq. (11) with the results of Monte Carlo simulation. Pinch point singularities are clearly seen in both analytic and numerical calculations. It is the same fluctuations of $\mathbf{B}_1$ and $\mathbf{B}_2$ which are responsible for the characteristic pinch-point structures in the (spin) structure factor measured by neutron scattering [cf. Fig. 2b].

At finite temperature, we anticipate that the spin liquid SL II will be perturbatively stable against terms such as $\delta\mathcal{H}_{\pm\pm}$ [Eq. (6)], which retain the point-group symmetry of the lattice but lift the $U(1)$ symmetry of the spins.

<figure><img src="image_5.png" /><figcaption>Figure 4. U(1)×U(1) gauge structure of the spin liquid SL II, as demonstrated by pinch-points in equal-time structure factors $S_{\mathbf{B}_{\mu}}^{\alpha\beta}(\mathbf{q})$ [Eq. (12)]. Left half of panel: structure factor $S_{\mathbf{B}_1}^{xx}(\mathbf{q})$ as calculated from the theory Eq. (11), with $\lambda = 1$. Right half of panel: structure factor $S_{\mathbf{B}_1}^{xx}(\mathbf{q})$ as calculated in classical Monte Carlo simulations of $\mathcal{H}_{QSI}$ [Eq. (1)]. The pinch point centered on $\mathbf{q} = (0,0,0)$, follows from the zero-divergence conditions on the fields $\mathbf{B}_{\mu}$ [Eq. (9)]. Simulations were carried out for a cubic cluster of $N = 8192$ spins, with $J^{\pm}/J_{zz} = -1$, $T = 0.01J_{zz}$, as described in Appendix E.</figcaption></figure>

In this case the free energy will be modified :

$$\mathcal{F}_{\text{SLII}} \rightarrow \mathcal{F}_{\text{SLII}} + \delta\mathcal{F}_{\text{SLII}} \quad (13)$$

$$\delta\mathcal{F}_{\text{SLII}} = \frac{T}{V} \int d^3\mathbf{r} \, \lambda' \left\{ (B_1^x)^2 - \frac{1}{2} [(B_1^y)^2 + (B_1^z)^2] - (B_2^x)^2 + \frac{1}{2} [(B_2^y)^2 + (B_2^z)^2] - \sqrt{3} [B_1^z B_2^z - B_1^y B_2^y] \right\} \quad (14)$$

This form of free energy will still lead to pinch points in $S_{\mathbf{B}_{\mu}}^{\alpha\beta}(\mathbf{q})$ and $S(\mathbf{q})$, but these will take on a more anisotropic character.

## IV. SPIN DYNAMICS IN THE SPIN-NEMATIC PHASE

For temperatures, $T < T_Q$ the spin-liquid SL II gives way to hidden spin-nematic order, labelled $Q_{\perp}$ in Fig. 1. As far as the dipole moments of spins are concerned, the spin-nematic phase is disordered, and neutron scattering experiments would reveal exactly the same algebraic correlations as in the spin liquid SL II. However the pinch points in $S(\mathbf{q})$ [cf Fig. 2b] hide a great wealth of interesting spin excitations.

To better understand the dynamics of the spin-nematic phase, we have calculated the dynamical structure factor $S(\mathbf{q}, \omega)$, within a semi-classical molecular-dynamics (MD) simulation, using the methods described in [69]. Relevant definitions are given in Appendix C. For $\omega/J_{zz} \lesssim 0.2$, $S(\mathbf{q}, \omega)$ presents a featureless, non-dispersing continuum [Fig. 5a]. Relics of dispersing ex-
<figure><img src="image_6.png" /><figcaption>Figure 5. Spin dynamics in the phase with hidden spin-nematic order ($Q_{\perp}$). (a) Dynamical structure factor for spin degrees of freedom, $S(q, \omega)$, showing gapless continuum of excitations at low energies. (b) Dynamical susceptibility for fluctuations of quadrupole moments, $\chi_{Q_{\perp}}(q, \omega)$. Inset: details of correlations at low energies near the zone center $q_{\parallel} = (2, 2, 2)$, showing a linearly-dispersing Goldstone mode at low energies. Results are taken from molecular-dynamics simulations of $H_{QSI}$ [Eq. (1)] for a cluster of $N = 65536$ spins, with $J_{\perp}/J_{zz} = -1.0$, $T/J_{zz} = 0.002$.</figcaption></figure>

citations are visible in $S(q, \omega)$ at higher energies, but these are explicitly not Goldstone modes, and have nothing to do with the hidden spin-nematic order. Results for $S(q, \omega)$ in the spin-nematic phase are very similar to those found in the spin liquid SL II for $T > T_Q$.

Incoherent, non-dispersing structure of the type shown in Fig. 5a is reminiscent of theoretical predictions [70–72] and experimental measurements [38, 73, 74], for a wide range of different spin liquids. The presence of a non-dispersing continuum reflects the fact that, unlike conventional spin waves (magnons), single elementary excitations of a spin liquid cannot be created by local processes. It follows that, when a neutron scatters from a spin liquid, the energy, momentum and angular momentum (spin) transferred are not absorbed by a single excitation with a well-defined energy and momentum, but rather shared between multiple excitations [75]. This scenario also holds for the spin liquid, SL II. However the fact that $S(q, \omega)$ only records dipolar spin correlations obscures another important fact — the spin-nematic order which is present for $T < T_Q$ which breaks a continuous, U(1), symmetry of the Hamiltonian. And, by Goldstone’s theorem, it must, therefore, also support gapless Goldstone modes.

In order to resolve this conundrum, it is necessary to examine the dynamical correlations of the quadrupole moments of spin. In Fig. 5b we present MD simulation results for the dynamical susceptibility $\chi_{Q_{\perp}}(q, \omega)$, which measures fluctuations of the order parameter $Q_{\perp}$ [Eq. (2)]. A sharp excitation, with dispersion

$$ \omega \approx v_Q |q - q_{\parallel}| \quad (15) $$

can now be resolved near to the zone centers with $q_{\parallel} = (0, 0, 0)$, $(1, 1, 1)$, $(2, 2, 2)$. These are the same zone centers for which the Bragg peaks associated with the hidden spin-nematic order $Q_{\perp}$ would occur in a quadrupolar structure factor, which might, in principle, be measured in resonant X-ray experiments [76].

It is interesting to compare these results with the behavior expected of a quantum spin nematic. Field-theoretic analysis [62, 77–79], based on the symmetry of the order parameter, predicts that spin-nematics support gapless Goldstone modes, visible in $\chi_{Q_{\perp}}(q, \omega)$. This Goldstone mode has dispersion $\omega \propto |q|$ [cf. Eq. (15)], and at zero temperature the associated intensity diverges as $\sim 1/\omega$ for $\omega \to 0$ [62]. In the present case, we find a linearly-dispersing Goldstone mode [Eq. (15)], with intensity which diverges for $\omega \to 0$ [cf inset to Fig. 5b].
However, the form of divergence observed in MD simu-
lations is ~1/ω². This follows from the fact that sim-
ulations are performed at finite temperature and probe
thermal rather than quantum fluctuations. Further de-
tails of the spin dynamics in the spin-nematic phase, and
specifically the characterization of the Goldstone mode
are given in Appendix F.

V. DISCUSSION

Spin–liquids and spin–nematics are both unconven-
tional phases expected to arise in certain frustrated mag-
nets, for which there has been a long experimental search.
This makes the possibility of realizing both of these in the
same system rather an enticing one. We therefore turn
to discuss the possibility of observing the SLII and Q⊥
phases in rare-earth pyrochlore magnets.

It is worth noting at this point that spin–nematic or-
der of the $Q_\perp$ phase is distinct from the multipolar orders
proposed as possible ground states of non-Kramers ions
($\text{Pr}^{3+}$ [8, 17, 38], $\text{Tb}^{3+}$ [8, 33, 39, 80]) and of dipole-
octupole Kramers doublets ($\text{Nd}^{3+}$ [44], $\text{Ce}^{3+}$ [81]). In
those cases the multipolar order is actually a “classical”
ordering of the transverse part of the pseudospins $\text{S}_i$.
This ordering must then be interpreted as a multipolar
order due to the symmetry properties of the crystal field
doublet which $\text{S}_i$ describes. Where multipole order of
this kind occurs, experiments which probe the dynamics
of dipoles will see a gapped response and a sharp exci-
tation spectrum. By contrast, in the case developed in
this Article, the dipole moments remain in an essentially
liquid like state with strong fluctuations at low temper-
ature and a broad, gapless response coexisting with the
hidden nematic order [Fig. 5a].

Where then might we observe this unusual magnetic
state? Further experimental work will be necessary to
definitively answer this question, but there are already
a few trails to follow. In particular, the Pr-based py-
rochlores have the recommended single-ion and inter-
action anisotropies [8, 17, 38]. Coupling parameters
of Pr₂Zr₂O₇ for example have been suggested to sit in
the AF⊥ phase of Fig. 1 [38], although it seems that
the coupling of structural disorder to the non-Kramers
doublets plays a significant role [29]. Since chemi-
cal pressure has already proven to be a useful tool to
move a family of compounds across a phase diagram
[36, 42, 67, 82, 83], Pr₂X₂O₇ (X=Sn,Hf,Pb) are promis-
ing candidates to investigate, with ferromagnetic corre-
lations consistent with positive *J*<sub>zz</sub> and no dipole order
yet observed [22, 27, 28, 35, 84].

The notion of hidden order also resonates with the elu-
sive physics of Yb-based pyrochlores. As far as we know,
Yb pyrochlores lie in a different regime of magnetic in-
teractions than the $\mathcal{H}_{\text{QSI}}$ model of Eq. (1) [10, 85, 86].
The properties of Yb₂Ti₂O₇ in particular seem to be
connected with a different phase boundary than the one
associated with SLII [42, 67]. That being said, the sim-

ilarities between our results and the Yb-pyrochlores are
striking: a gapless continuum of spin excitations, obliv-
ious to the transition [36, 87, 88] [Fig. 5(a)], and robust
in temperature up to a broad feature in specific heat [36]
(here between SL III and SL II). And while the magnetic
order in Yb-pyrochlores is, at least partially, an order
of dipolar moments [23, 31, 32, 89–91], recent experi-
ments have indicated that the primary order parameter
may be “hidden”, and distinct from a standard dipole
order [36]. Thus, although the specific case developed
in this manuscript probably does not apply to the Yb-
pyrochlores, related physics may be at play.

Furthermore, since nematic order does not prevent
dipole disorder, our work provides a model Hamilto-
nian for the co-existence between emergent gauge fields
and long-range order. Hence, the $Q_\perp$ phase offers an
unconventional counterpart to the Coulombic ferromag-
net [16, 92], and spin-ice fragmentation [93] that has re-
cently been observed in Nd₂Zr₂O₇ [37, 94] and Ho₂Ir₂O₇
[95] materials.

We also note that many other magnetic systems out-
side the rare earth oxides R₂X₂O₇ feature moments lo-
cated on a pyrochlore lattice. Of particular interest are
materials such as NaCaCo₂F₇ and NaSrCo₂F₇ [96, 97]
which boast XY like interactions with much higher en-
ergy scales than observed in the rare-earth oxides. If such
a case could be found with frustrated transverse coupling
$J_± < 0$ then it would render the physics discussed here
accessible at a much more amenable temperature range.

VI. SUMMARY AND CONCLUSIONS

"Quantum spin ice", in which magnetic ions on a pyrochlore lattice interact through highly-anisotropic exchange interactions, have become an important paradigm in the search for quantum spin liquids. In this Article we have used large-scale Monte-Carlo simulation to explore the physics of the minimal model of a quantum spin ice $\mathcal{H}_{\text{QSI}}$ [1]. We find that this model has far more to offer than spin ice alone, supporting three distinct types of spin liquid, each with a different emergent gauge symmetry. These include a completely new form of spin liquid, described by a $U(1) \times U(1)$ gauge theory, which exhibits a hidden spin-nematic order at low temperatures [Fig. 1]. So far as experiment is concerned, the main lesson of these results is that "quantum spin-ice" materials can play host to a great many different spin-liquid and (hidden-)order phases — crudely-put : pinch-points do not imply spin ice [Fig. 2]. The existence of a sharp Goldstone mode in the nematic phase $Q_\perp$ also serves as a salutary reminder that broad, non-dispersing continua of excitations can hide a multitude of secrets [Fig. 5].

From the theoretical point of view, this work opens an
interesting new perspective on quantum spin liquids on
the pyrochlore lattice. The effect of quantum fluctua-
tions on the novel ordered and spin-liquid phases shown
in Fig. 1 for $J_\pm < 0$ remains an open question. Experi-
ence with QMC simulation of $\mathcal{H}_{QSI}$ [Eq. (1)] for $J_{\pm} > 0$ suggests that quantitative values of the crossover temperature $T_2^*$ and $T_3^*$ may be substantially renormalized, but that the qualitative structure of the phase diagram should remain the same down to very low temperatures [13, 20, 21]. The high-symmetry point, $J_{\pm}/J_{zz} = -1/2$ is also a high-symmetry point for quantum spins, and so remains the anchor for the spin liquid SL III (pHAF) None the less, the fate of this $U(1) \times U(1) \times U(1)$ spin liquid for quantum spins at $T = 0$ remains an open question [50, 51, 53, 98]. And, to the best of our knowledge, quantum analogues of the new spin liquid, SL II, which has a $U(1) \times U(1)$ gauge structure, remain unexplored. However it seems reasonable to speculate that quantum effects will enhance, rather than suppress, the fluctuations which drive SL II and SL III, and that the phase $Q_{\perp}$ will survive as hidden quantum spin-nematic order, within a quantum spin liquid. All of these questions open exciting avenues for future research.

## ACKNOWLEDGEMENTS

We are pleased to acknowledge stimulating discussions with Karim Essafi. This work was supported by the Theory of Quantum Matter Unit of the Okinawa Institute of Science and Technology Graduate University (OIST). Numerical calculations were carried out using HPC facilities provided by OIST.

## Appendix A: Definition of local coordinate frame

We describe the local-coordinate frame which is defined for four spins on a pyrochlore tetrahedron $S_0, S_1, S_2, S_3$ occupying positions

$$ \begin{aligned} \mathbf{r}_0 &= \frac{a}{8} (1, 1, 1) & \mathbf{r}_1 &= \frac{a}{8} (1, -1, -1) \\ \mathbf{r}_2 &= \frac{a}{8} (-1, 1, -1) & \mathbf{r}_3 &= \frac{a}{8} (-1, -1, 1), \end{aligned} \quad (A1) $$

where $a$ is the length of a cubic, 16-site unit cell of the pyrochlore lattice.

The pseudospins in the global, crystal, coordinate frame $S_i$ relate to the pseudospins in the local frame $S_i$ [Eq. (1)] as

$$ \mathbf{S}_i = \mathbf{x}_i^{\text{local}} \mathbf{S}_i^x + \mathbf{y}_i^{\text{local}} \mathbf{S}_i^y + \mathbf{z}_i^{\text{local}} \mathbf{S}_i^z \quad (A2) $$

Where

$$ \begin{aligned} \mathbf{z}_0^{\text{local}} &= \frac{1}{\sqrt{3}}(1, 1, 1) & \mathbf{z}_1^{\text{local}} &= \frac{1}{\sqrt{3}}(1, -1, -1) \\ \mathbf{z}_2^{\text{local}} &= \frac{1}{\sqrt{3}}(-1, 1, -1) & \mathbf{z}_3^{\text{local}} &= \frac{1}{\sqrt{3}}(-1, -1, 1), \end{aligned} \quad (A3) $$

$$ \begin{aligned} \mathbf{x}_0^{\text{local}} &= \frac{1}{\sqrt{6}}(-2, 1, 1) & \mathbf{x}_1^{\text{local}} &= \frac{1}{\sqrt{6}}(-2, -1, -1) \\ \mathbf{x}_2^{\text{local}} &= \frac{1}{\sqrt{6}}(2, 1, -1) & \mathbf{x}_3^{\text{local}} &= \frac{1}{\sqrt{6}}(2, -1, 1), \end{aligned} \quad (A4) $$

and

$$ \begin{aligned} \mathbf{y}_0^{\text{local}} &= \frac{1}{\sqrt{2}}(0, -1, 1) & \mathbf{y}_1^{\text{local}} &= \frac{1}{\sqrt{2}}(0, 1, -1) \\ \mathbf{y}_2^{\text{local}} &= \frac{1}{\sqrt{2}}(0, -1, -1) & \mathbf{y}_3^{\text{local}} &= \frac{1}{\sqrt{2}}(0, 1, 1). \end{aligned} \quad (A5) $$

## Appendix B: Details of the numerical determination of the phase diagram

The phase diagram shown in Fig. 1 was extracted from classical Monte Carlo (MC) simulations of the quantum spin ice model, $\mathcal{H}_{QSI}$ [Eq. (1)]. Spins were treated as classical vectors with fixed length $|S_i| = \frac{1}{2}$. These simulations were carried out for a cubic cluster of 8192 spins, using a single spin flip algorithm combined with simulated annealing, parallel tempering and over-relaxation. The phase diagram is obtained using 75000 simulated annealing steps using a stepwise decrease of temperature starting from $T = 10J_{zz}$ down to the target temperature. Each annealing step consists of 10 Monte Carlo steps (a Monte Carlo step consists of a full sweep of the lattice combined with over-relaxation). The simulated annealing is followed by 1000 parallel tempering steps with

<figure><img src="image_7.png" /><figcaption>Figure 6. (Color online). Crossover in the magnetic susceptibility $\chi(T)$, from a traditional Curie law in the high-temperature paramagnet, to a low-temperature Curie law in the spin liquids, as seen by different plateaux in the function $T\chi$, plotted as a function of $\log(T)$. The crossover temperature $T_3^*/J_{zz} \approx 0.3$ (red dot) is estimated from the point of inflection of $T\chi$. The extraction of the crossover temperature $T_2^*$ is explained in Fig. 7. Results are taken from classical Monte Carlo simulations of $\mathcal{H}_{QSI}$ [Eq (1)], for a cubic cluster of $N = 8192$ spins, with $J_{\pm}/J_{zz} = -1$.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 7. Thermodynamics of the QSI in the region of spin-nematic order. (a) Specific heat $c_V(T)$, showing an upturn followed by a shallow maximum in the region of the crossover into the spin liquid SL II at $T_2^*/J_{zz} \approx 10^{-1}$, and small peak associated with the onset of spin-nematic order at $T_Q/J_{zz} \approx 10^{-2}$. (b) Correlation function $T_{\chi_{T_1\text{ice}}}(T)$ used to determine the crossover temperature into the spin-liquid SLII. $\chi_{T_1\text{ice}}$ is the susceptibility of the field $m_{T_1\text{ice}}$ defined in Appendix D. The crossover temperature $T_2^*$ is defined by the point at which the quantity $T_{\chi_{T_1\text{ice}}}(T)$ drops below its infinite temperature limit [Eq. B1] (b) Norm of the order parameter $Q_{\perp}(T)$ [cf. Eq. (2)], showing continuous phase transition into phase with hidden spin-nematic order at $T_Q/J_{zz} \approx 10^{-2}$. (c) Order-parameter susceptibility $\chi_{Q_{\perp}}(T)$ [Eq. (B2)], used to estimate the transition temperature $T_Q$. The crossover temperature $T_3^*$ has been extracted as explained in Fig. 6. Results are taken from classical Monte Carlo simulation of $\mathcal{H}_{\text{QSI}}$ [Eq. (1)], for cubic clusters of $N = 1024$, $N = 8192$ and $N = 27648$ spins, with $J_{\pm}/J_{zz} = -1$.</figcaption></figure>

500 Monte Carlo steps in between, and then by 200000 Monte Carlo steps for thermalization at fixed temperature. Measurements consist of 200000 samples separated by 10 Monte Carlo steps and combined to parallel tempering every 50 measures. We use 256 different replicas with temperature set in linear scale for $J_{\pm}/J_{zz} > -\frac{1}{2}$ and 256 temperatures in logarithmic scale for $J_{\pm}/J_{zz} \le -\frac{1}{2}$.

The phase boundary of the antiferromagnetically ordered ($AF_{\perp}$) phase, $T_N$, was extracted from the susceptibility of the relevant order parameter, $m_E$, as defined in Appendix D.

The crossover scale for the spin-ice regime (SL I), $T_1^*$, was estimated from the Schottky-like peak in the heat capacity.

The crossover scale $T_3^*$ for the spin-liquid SLIII was estimated from the Curie-Law crossover shown in Fig. 6.

For $J_{\pm} < -\frac{1}{2}$, the crossover scale $T_2^*$ is associated with a weakening of the correlations of the local z-components of the spins. This can be observed by measuring the susceptibility, $\chi_{T_1\text{ice}}(T)$, of the field $m_{T_1\text{ice}}$, defined in Appendix D]. Decreasing the temperature for $-1 < \frac{J_{\pm}}{J_{zz}} < -0.5$ the quantity $T_{\chi_{T_1\text{ice}}}(T)$ first increases during the crossover from the paramagnet to SLIII and then drops as the system enters SLII. We define the crossover temperature $T_2^*$ as the point at which the quantity $T_{\chi_{T_1\text{ice}}}(T)$ drops below its infinite temperature value

$$T_2^* \chi_{T_1\text{ice}}(T_2^*) = \lim_{T \to \infty} T \chi_{T_1\text{ice}}(T) \quad (B1)$$

This is illustrated in Fig. 7b.

The spin-nematic ordering temperature $T_Q$ is estimated from the peak in the order-parameter susceptibility

$$\chi_{Q_{\perp}} = \frac{N}{T} (\langle Q_{\perp}^2 \rangle - \langle |Q_{\perp}| \rangle^2). \quad (B2)$$

[cf. Fig. 7d], where the relevant order parameter $Q_{\perp}$ [cf. Fig. 7c], is defined through Eq. (2).

Fig. 7 is obtained using 300 temperatures in logarithmic scale covering 3 orders of magnitude, parallel tempering every 100 Monte Carlo steps, simulated annealing and thermalization at temperature $T$ for 100000 Monte Carlo steps each. Measurements consist of 100000 different samples with 10 Monte Carlo steps between each sample. Error bars were estimated by comparing the results of three independent runs of the simulation.

## Appendix C: Definitions of dynamical structure factors

In Fig. 2 we show predictions for neutron scattering experiments, based on the equal-time (i.e. energy-integrated) structure factor

$$S(\mathbf{q}) = \int d\omega S(\mathbf{q}, \omega), \quad (C1)$$
where the dynamical structure factor $S(\mathbf{q}, \omega)$ is defined through

$$S(\mathbf{q}, \omega) = \sum_{\alpha\beta} \left( \delta_{\alpha\beta} - \frac{q_{\alpha}q_{\beta}}{q^2} \right) \langle m^{\alpha}(-\mathbf{q}, \omega) m^{\beta}(\mathbf{q}, \omega) \rangle$$

$$m^{\alpha}(\mathbf{q}, \omega) = \sum_{i\beta\gamma} R_{i}^{\alpha\beta} g_{i}^{\beta\gamma} \left( \int S_{i}^{\gamma}(t) e^{i\omega t} dt \right) e^{i\mathbf{q} \cdot \mathbf{r}_{i}}, \quad (C2)$$

and the $g_i^{\beta\gamma}$ is the g-tensor written in the local coordinate frame [42]. For simplicity, we have here taken $g_i^{\beta\gamma} = 2\delta_{\beta\gamma}$ for all of the calculations in this paper. $R_i^{\alpha\beta}$ is a rotation matrix which rotates from the local coordinate frame on site $i$, to the global, crystal coordinate frame. The definition of the local coordinate frame is given in Appendix A. Results for $S(\mathbf{q})$ are shown in the left half-panels of Fig. 2. These results were taken from classical MC simulations of $\mathcal{H}_{QSI}$ at a given temperature, with further averaging provided by numerically integrating the semi-classical equations of motion for the spins. This secondary molecular-dynamics (MD) simulation was carried out using methods described in Ref. [69].

It is also useful to decompose the structure factor into the spin-flip (SF) and non spin-flip (NSF) channels measured in polarised neutron-scattering experiments.

$$S_{\text{SF}}(\mathbf{q}) = \frac{1}{q^2} \int d\omega \langle |\mathbf{m}(\mathbf{q}, \omega) \cdot (\hat{\mathbf{n}} \times \mathbf{q})|^2 \rangle \\ S_{\text{NSF}}(\mathbf{q}) = \int d\omega \langle |\mathbf{m}(\mathbf{q}, \omega) \cdot \hat{\mathbf{n}}|^2 \rangle, \quad (C3)$$

where $\hat{\mathbf{n}}$ is the direction of polarization of the neutron magnetic moment. Following Fennell et al. [54], we take $\hat{\mathbf{n}} = (1, -1, 0)/\sqrt{2}$. Simulation results for $S_{\text{SF}}(\mathbf{q})$ and $S_{\text{NSF}}(\mathbf{q})$ are shown in the right half-panels of Fig. 2.

We have also used MD simulation to calculate the dynamical structure factor $S(\mathbf{q}, \omega)$. Results for $S(\mathbf{q}, \omega)$ within the spin-nematic phase of the quantum spin ice model are shown in Fig. 5a. Further details of the calculation of dynamical properties can be found in Appendix F.

## Appendix D: Definitions of local order parameter fields

The definitions of the local order parameter fields $m_\lambda$ which appear in the theory of the spin liquid SLII [Section III] are given in Table I.

Here we give the definitions in terms of the spins written in the local coordinate frame $S_i$ (defined in Appendix A), cf. Ref. 42 where definitions are given in the global, crystal basis.

<table><caption>Table I. Order-parameter fields $m_\lambda$, derived from irreducible representations (irreps) of the tetrahedral point-group $T_d$. Spin components $S_i = (S_i^x, S_i^y, S_i^z)$ are written in the local frame of the magnetic ions, see Appendix A for a definition of this coordinate frame. The convention for the labelling of the spins with an tetrahedron is given in Appendix A.</caption><thead><tr><th></th><th>Definition in terms of spins within tetrahedron</th></tr></thead><tbody><tr><td>$m_{A_2}$</td><td>$\frac{1}{2}(S_0^x + S_1^x + S_2^x + S_3^x)$</td></tr><tr><td>$m_E$</td><td>$\frac{1}{2}(S_0^x + S_1^x + S_2^x + S_3^x)$</td></tr><tr><td>$m_{T_1,ice}$</td><td>$\frac{1}{2}(S_0^z - S_1^z - S_2^z - S_3^z)$</td></tr><tr><td>$m_{T_1,planar}$</td><td>$\frac{1}{4}(-S_0^x + \sqrt{3}S_0^y + S_1^x - \sqrt{3}S_1^y - S_2^x + \sqrt{3}S_2^y + S_3^x - \sqrt{3}S_3^y)$</td></tr><tr><td>$m_{T_2}$</td><td>$\frac{1}{4}(-S_0^x - \sqrt{3}S_0^y + S_1^x - \sqrt{3}S_1^y + S_2^x + \sqrt{3}S_2^y - S_3^x - \sqrt{3}S_3^y)$</td></tr></tbody></table>

## Appendix E: Numerical simulation of the correlations of the flux

Values of the flux field $\mathbf{B}_\mu(\mathbf{r})$ are calculated for each tetrahedron $\mathbf{r}$ according to Eq. (8) and the definitions of $m_\lambda$ given in Table I.

The tetrahedra of the pyrochlore lattice may be divided into two sets A and B. The centres of each set of tetrahedra each form an FCC lattice.

To calculate $\mathbf{S}_{\mathbf{B}_\mu}^{\alpha\beta}$ we use Eq. (12) where $\mathbf{B}_\mu(\mathbf{q})$ is defined as the lattice Fourier transform of $\mathbf{B}_\mu(\mathbf{r})$ over only the A sublattice of tetrahedra.

$$\mathbf{B}_{\mu}(\mathbf{q}) = \sqrt{\frac{1}{N_{\text{uc}}} \sum_{\mathbf{r} \in \mathbf{r}_{\text{A}}} \exp(-i\mathbf{q} \cdot \mathbf{r}_{\text{A}}) \mathbf{B}_{\mu}(\mathbf{r})} \quad (\text{E1})$$

where $N_{\text{uc}}$ is the number of unit cells in the system.

Simulations were carried out using local spin updates, augmented by over-relaxation, within a parallel tempering scheme with 300 temperatures distributed on a log scale between $T = 0.003 J_{zz}$ and $T = 0.1 J_{zz}$. Thermalisation was accomplished through a process of simulated annealing, with $10^4$ Monte Carlo steps (MCs) of annealing from high temperature to temperature $T$, followed by $10^4$ MCs of thermalization at temperature $T$, and $10^5$ MCs of measurements at temperature $T$. Spin configurations were sampled every 100 MCs during the measurements, giving an ensemble of 1000 samples.

## Appendix F: Dynamics of excitations in the spin-nematic phase

To study the Goldstone mode associated to the development of spin-nematic order, we calculate the dynamical correlation function.

$$\chi_{Q\perp}(\mathbf{q}, \omega) = \langle |\delta\mathbf{Q}_{\perp}(\mathbf{q}, \omega)|^2 \rangle, \quad (\text{F1})$$
<figure><img src="image_9.png" /><figcaption>Figure 8. Dispersion and intensity of the Goldstone mode in the phase with hidden spin-nematic order, as shown in Fig. 5b. (a) Dispersion $\epsilon(q)$ of low-energy peak in $\chi_{Q\perp}(\mathbf{q}, \omega)$ showing the expected behaviour $\epsilon(q) = q$ at small $q$. (b) Intensity $I(q)$ of the peak as a function of momentum $q$. The dashed line shows the expected behaviour at finite temperature, $I(q) \propto 1/q^2$. Results are taken from molecular-dynamics simulations of a cluster of $N = 65536$ spins, for $J_{\pm}/J_{zz} = -1.0$, $T = 0.002J_{zz}$. Momentum $q$ is measured relative to $\mathbf{q} = (0, 0, 0)$.</figcaption></figure>

where fluctuations of spin-nematic order are given by

$$ \delta \mathbf{Q}_{\perp}(\mathbf{q}, \omega) = \sum_{i} \int dt \left[ \mathbf{Q}_{\perp}(\mathbf{r}_i, t) - \bar{\mathbf{Q}}_{\perp}(t) \right] e^{i\omega t} e^{i\mathbf{q} \cdot \mathbf{r}_i} \\ \bar{\mathbf{Q}}_{\perp}(t) = \frac{1}{N} \sum_{i} \mathbf{Q}_{\perp}(\mathbf{r}_i, t), \qquad (\text{F2}) $$

and the order parameter $\mathbf{Q}_{\perp}(\mathbf{r}_i, t)$ is defined through Eq. (2).

$\chi_{Q\perp}(q, \omega)$ is calculated numerically from 200 sample configurations extracted from Monte Carlo simulations on a system of linear size $L = 16$. We used 20000 steps for the simulated annealing spaced by 10 Monte Carlo steps between each simulated annealing step. The other parameters for the thermalization and parallel tempering are identical to the parameters used to calculate the phase diagram [Appendix B].

The ensemble of configurations obtained from Monte Carlo is then evolved in time according to the equation of motion,

$$ \frac{dS_i}{dt} = H_i(t) \times S_i(t) \qquad (\text{F3}) $$

where

$$ H_i(t) = \sum_{j \in nn} J_{ij} \cdot S_j(t) \qquad (\text{F4}) $$

is the effective exchange field acting on site $i$, $J_{ij}$ is the anisotropic exchange interaction tensor and the sum in Eq. (F4) runs over the neighbors of $i$. The numerical integration of this nonlinear equation of motion proceeds as described in Ref. [69].

In Fig 8 we plot the dispersion of the Goldstone mode found in molecular dynamics (MD) simulations of $H_{QSI}$ within the spin-nematic phase for $J_{\pm}/J_{zz} = -1$. The dispersion was extracted from the position of low-energy dispersing peak in $\chi_{Q\perp}(\mathbf{q}, \omega)$, as shown in the inset to Fig. 5b.

[1] Patrick A. Lee, “An end to the drought of quantum spin liquids,” Science 321, 1306–1307 (2008).

[2] Leon Balents, “Spin liquids in frustrated magnets,” Nature 464, 199 (2010).

[3] Jason S. Gardner, Michel J. P. Gingras, and John E. Greedan, “Magnetic pyrochlore oxides,” Reviews of Modern Physics 82, 53–107 (2010).

[4] Steven T. Bramwell and Michel J. P. Gingras, “Spin Ice State in Frustrated Magnetic Pyrochlore Materials,” Science 294 (2001).

[5] Claudio Castelnovo, Roderich Moessner, and S.L. Sondhi, “Spin Ice, Fractionalization, and Topological Order,” Annu. Rev. Condens. Matter Phys. 3, 35 (2012).

[6] S. H. Curnoe, “Quantum spin configurations in $Tb_2Ti_2O_7$,” Phys. Rev. B 75, 212404 (2007).

[7] Hamid R. Molavian, Michel J. P. Gingras, and Benjamin Canals, “Dynamically Induced Frustration as a Route to a Quantum Spin Ice State in $Tb_2Ti_2O_7$ via Virtual Crys-

tal Field Excitations and Quantum Many-Body Effects,” Phys. Rev. Lett. 98, 157204 (2007).

[8] Shigeki Onoda and Yoichi Tanaka, “Quantum melting of spin ice: Emergent cooperative quadrupole and chirality,” Phys. Rev. Lett. 105, 047201 (2010).

[9] Shigeki Onoda and Yoichi Tanaka, “Quantum fluctuations in the effective pseudospin-$\frac{1}{2}$ model for magnetic pyrochlore oxides,” Phys. Rev. B 83, 094411 (2011).

[10] Kate A. Ross, Lucile Savary, Bruce D. Gaulin, and Leon Balents, “Quantum excitations in quantum spin ice,” Phys. Rev. X 1, 021002 (2011).

[11] M J P Gingras and P A McClarty, “Quantum spin ice: a search for gapless quantum spin liquids in pyrochlore magnets,” Reports on Progress in Physics 77, 056501 (2014).

[12] Michael Hermele, Matthew P. A. Fisher, and Leon Balents, “Pyrochlore photons: The U(1) spin liquid in a $S=\frac{1}{2}$ three-dimensional frustrated magnet,” Phys. Rev.
B 69, 064404 (2004).

[13] Argha Banerjee, Sergei V. Isakov, Kedar Damle, and Yong Baek Kim, “Unusual liquid state of hard-core bosons on the pyrochlore lattice,” Phys. Rev. Lett. 100, 047208 (2008).

[14] Nic Shannon, Olga Sikora, Frank Pollmann, Karlo Penc, and Peter Fulde, “Quantum ice: A quantum monte carlo study,” Phys. Rev. Lett. 108, 067204 (2012).

[15] Owen Benton, Olga Sikora, and Nic Shannon, “Seeing the light: Experimental signatures of emergent electromagnetism in a quantum spin ice,” Phys. Rev. B 86, 075154 (2012).

[16] Lucile Savary and Leon Balents, “Coulombic quantum liquids in spin-1/2 pyrochlores,” Phys. Rev. Lett. 108, 037202 (2012).

[17] SungBin Lee, Shigeki Onoda, and Leon Balents, “Generic quantum spin ice,” Phys. Rev. B 86, 104412 (2012).

[18] Zhihao Hao, Alexandre G. R. Day, and Michel J. P. Gingras, “Bosonic many-body theory of quantum spin ice,” Phys. Rev. B 90, 214430 (2014).

[19] P. A. McClarty, O. Sikora, R. Moessner, K. Penc, F. Pollmann, and N. Shannon, “Chain-based order and quantum spin liquids in dipolar spin ice,” Phys. Rev. B 92, 094418 (2015).

[20] Yasuyuki Kato and Shigeki Onoda, “Numerical evidence of quantum melting of spin ice: Quantum-to-classical crossover,” Phys. Rev. Lett. 115, 077202 (2015).

[21] Nic Shannon, “Spin ice,” (Springer, 2017) Chap. “Quantum Monte Carlo simulations of quantum spin ice”.

[22] H. D. Zhou, C. R. Wiebe, J. A. Janik, L. Balicas, Y. J. Yo, Y. Qiu, J. R. D. Copley, and J. S. Gardner, “Dynamical Spin Ice: Pr₂Sn₂O₇,” Phys. Rev. Lett. 101, 227204 (2008).

[23] Lieh-Jeng Chang, Shigeki Onoda, Yixi Su, Ying-Jer Kao, Ku-Ding Tsuei, Yukio Yasui, Kazuhisa Kakurai, and Martin Richard Lees, “Higgs transition from a magnetic Coulomb liquid to a ferromagnet in Yb₂Ti₂O₇,” Nat. Commun. 3, 992 (2012).

[24] T. Fennell, M. Kenzelmann, B. Roessli, M. K. Haas, and R. J. Cava, “Power-Law Spin Correlations in the Pyrochlore Antiferromagnet Tb₂Ti₂O₇,” Phys. Rev. Lett. 109, 017201 (2012).

[25] K. Kimura, S. Nakatsuji, J.-J. Wen, C. Broholm, M. B. Stone, E. Nishibori, and H. Sawa, “Quantum fluctuations in spin-ice-like Pr₂Zr₂O₇,” Nat Commun 4 (2013).

[26] Romain Sibille, Elsa Lhotel, Vladimir Pomjakushin, Chris Baines, Tom Fennell, and Michel Kenzelmann, “Candidate Quantum Spin Liquid in the Ce³⁺ Pyrochlore Stannate Ce₂Sn₂O₇,” Phys. Rev. Lett. 115, 097202 (2015).

[27] Romain Sibille, Elsa Lhotel, Monica Ciomaga Hatnean, Geetha Balakrishnan, Björn Fåk, Nicolas Gauthier, Tom Fennell, and Michel Kenzelmann, “Candidate quantum spin ice in the pyrochlore Pr₂Hf₂O₇,” Phys. Rev. B 94, 024436 (2016).

[28] V. K. Anand, L. Opherden, J. Xu, D. T. Adroja, A. T. M. N. Islam, T. Herrmannsdörfer, J. Hornung, R. Schönemann, M. Uhlarz, H. C. Walker, N. Casati, and B. Lake, “Physical properties of the candidate quantum spin-ice system Pr₂Hf₂O₇,” Phys. Rev. B 94, 144415 (2016).

[29] J.-J. Wen, S. M. Koohpayeh, K. A. Ross, B. A. Trump, T. M. McQueen, K. Kimura, S. Nakatsuji, Y. Qiu, D. M. Pajerowski, J. R. D. Copley, and C. L. Broholm, “Disordered Route to the Coulomb Quantum Spin Liquid: Random Transverse Fields on Spin Ice in Pr₂Zr₂O₇,” Phys. Rev. Lett. 118, 107206 (2017).

[30] P. Dalmas de Réotier, A. Yaouanc, L. Keller, A. Cervellino, B. Roessli, C. Baines, A. Forget, C. Vaju, P. C. M. Gubbens, A. Amato, and P. J. C. King, “Spin Dynamics and Magnetic Order in Magnetically Frustrated Tb₂Sn₂O₇,” Phys. Rev. Lett. 96, 127202 (2006).

[31] Z. L. Dun, E. S. Choi, H. D. Zhou, A. M. Hallas, H. J. Silverstein, Y. Qiu, J. R. D. Copley, J. S. Gardner, and C. R. Wiebe, “Yb₂Sn₂O₇: A magnetic Coulomb liquid at a quantum critical point,” Phys. Rev. B 87, 134408 (2013).

[32] A. Yaouanc, P. Dalmas de Réotier, P. Bonville, J. A. Hodges, V. Glazkov, L. Keller, V. Sikolenko, M. Bartkowiak, A. Amato, C. Baines, P. J. C. King, P. C. M. Gubbens, and A. Forget, “Dynamical Splayed Ferromagnetic Ground State in the Quantum Spin Ice Yb₂Sn₂O₇,” Phys. Rev. Lett. 110, 127207 (2013).

[33] T. Taniguchi, H. Kadowaki, H. Takatsu, B. Fåk, J. Ollivier, T. Yamazaki, T. J. Sato, H. Yoshizawa, Y. Shimura, T. Sakakibara, T. Hong, K. Goto, L. R. Yaraskavitch, and J. B. Kycia, “Long-range order and spin-liquid states of polycrystalline Tb₂+xTi₂-xO₇+y,” Phys. Rev. B 87, 060408 (2013).

[34] A. M. Hallas, J. G. Cheng, A. M. Arevalo-Lopez, H. J. Silverstein, Y. Su, P. M. Sarte, H. D. Zhou, E. S. Choi, J. P. Attfield, G. M. Luke, and C. R. Wiebe, “Incipient Ferromagnetism in Tb₂Ge₂O₇: Application of Chemical Pressure to the Enigmatic Spin-Liquid Compound Tb₂Ti₂O₇,” Phys. Rev. Lett. 113, 267205 (2014).

[35] A. M. Hallas, A. M. Arevalo-Lopez, A. Z. Sharma, T. Munsie, J. P. Attfield, C. R. Wiebe, and G. M. Luke, “Magnetic frustration in lead pyrochlores,” Phys. Rev. B 91, 104417 (2015).

[36] A. M. Hallas, J. Gaudet, N. P. Butch, M. Tachibana, R. S. Freitas, G. M. Luke, C. R. Wiebe, and B. D. Gaulin, “Universal dynamic magnetism in Yb pyrochlores with disparate ground states,” Phys. Rev. B 93, 100403 (2016).

[37] S. Petit, E. Lhotel, B. Canals, M. Ciomaga Hatnean, J. Ollivier, H. Mutka, E. Ressouche, A. R. Wildes, M. R. Lees, and G. Balakrishnan, “Observation of magnetic fragmentation in spin ice,” Nat Phys 12, 746–750 (2016).

[38] S. Petit, E. Lhotel, S. Guitteny, O. Florea, J. Robert, P. Bonville, I. Mirebeau, J. Ollivier, H. Mutka, E. Ressouche, C. Decorse, M. Ciomaga Hatnean, and G. Balakrishnan, “Antiferroquadrupolar correlations in the quantum spin ice candidate Pr₂Zr₂O₇,” Phys. Rev. B 94, 165153 (2016).

[39] H. Takatsu, S. Onoda, S. Kittaka, A. Kasahara, Y. Kono, T. Sakakibara, Y. Kato, B. Fåk, J. Ollivier, J. W. Lynn, T. Taniguchi, M. Wakita, and H. Kadowaki, “Quadrupole Order in the Frustrated Pyrochlore Tb₂+xTi₂-xO₇+y,” Phys. Rev. Lett. 116, 217201 (2016).

[40] Shigeki Onoda, “Effective quantum pseudospin-1/2 model for Yb pyrochlore oxides,” Journal of Physics: Conference Series 320, 012065 (2011).

[41] Anson W. C. Wong, Zhihao Hao, and Michel J. P. Gingras, “Ground state phase diagram of generic xy pyrochlore magnets with quantum fluctuations,” Phys. Rev. B 88, 144402 (2013).
[42] Han Yan, Owen Benton, Ludovic Jaubert, and Nic Shannon, “Theory of multiple-phase competition in pyrochlore magnets with anisotropic exchange with application to Yb₂Ti₂O₇, Er₂Ti₂O₇, and Er₂Sn₂O₇,” Phys. Rev. B 95, 094422 (2017).

[43] Owen Benton, *Classical and quantum spin liquids on the pyrochlore lattice*, Ph.D. thesis, University of Bristol (2015).

[44] Yi-Ping Huang, Gang Chen, and Michael Hermele, “Quantum spin ices and topological phases from dipolar-octupolar doublets on the pyrochlore lattice,” Phys. Rev. Lett. 112, 167203 (2014).

[45] G. Chen, “The Spectral Periodicity of Spinon Continuum in Quantum Spin Ice,” arXiv:1704.02734.

[46] R. Moessner and J. T. Chalker, “Properties of a classical spin liquid: The heisenberg pyrochlore antiferromagnet,” Phys. Rev. Lett. 80, 2929–2932 (1998).

[47] R. Moessner and J. T. Chalker, “Low-temperature properties of classical geometrically frustrated antiferromagnets,” Phys. Rev. B 58, 12049–12062 (1998).

[48] S. V. Isakov, K. Gregor, R. Moessner, and S. L. Sondhi, “Dipolar spin correlations in classical pyrochlore magnets,” Phys. Rev. Lett. 93, 167204 (2004).

[49] C. L. Henley, “Power-law spin correlations in pyrochlore antiferromagnets,” Phys. Rev. B 71, 014424 (2005).

[50] Benjamin Canals and Claudine Lacroix, “Pyrochlore antiferromagnet: A three-dimensional quantum spin liquid,” Phys. Rev. Lett. 80, 2933 (1998).

[51] Benjamin Canals and Claudine Lacroix, “Quantum spin liquid: The heisenberg antiferromagnet on the three-dimensional pyrochlore lattice,” Phys. Rev. B 61, 1149 (2000).

[52] F. J. Burnell, Shoibal Chakravarty, and S. L. Sondhi, “Monopole flux state on the pyrochlore lattice,” Phys. Rev. B 79, 144432 (2009).

[53] Yuan Huang, Kun Chen, Youjin Deng, Nikolay Prokof'ev, and Boris Svistunov, “Spin-ice state of the quantum heisenberg antiferromagnet on the pyrochlore lattice,” Phys. Rev. Lett. 116, 177203 (2016).

[54] T. Fennell, P. P. Deen, A. R. Wildes, K. Schmalzl, D. Prabhakaran, A. T. Boothroyd, R. J. Aldus, D. F. McMorrow, and S. T. Bramwell, “Magnetic Coulomb Phase in the Spin Ice Ho₂Ti₂O₇,” Science 326, 415–417 (2009).

[55] Owen Benton, unpublished.

[56] Mathieu Taillefumier, Owen Benton, and Nic Shannon, In preparation.

[57] L. D. C. Jaubert, M. J. Harris, T. Fennell, R. G. Melko, S. T. Bramwell, and P. C. W. Holdsworth, “Topological-sector fluctuations and curie-law crossover in spin ice,” Phys. Rev. X 3, 011014 (2013).

[58] P. H. Conlon and J. T. Chalker, “Spin dynamics in pyrochlore heisenberg antiferromagnets,” Phys. Rev. Lett. 102, 237206 (2009).

[59] Henley, C. L., “The “Coulomb Phase” in Frustrated Systems,” Annu. Rev. Condens. Matter Phys. 1, 179–2010 (2010).

[60] Nic Shannon, Tsutomu Momoi, and Philippe Sindzinger, “Nematic order in square lattice frustrated ferromagnets,” Phys. Rev. Lett. 96, 027213 (2006).

[61] Nic Shannon, Karlo Penc, and Yukitoshi Motome, “Nematic, vector-multipole, and plateau-liquid states in the classical $O(3)$ pyrochlore antiferromagnet with bi-quadratic interactions in applied magnetic field,” Phys.

Rev. B 81, 184409 (2010).

[62] Andrew Smerald, Hiroaki T. Ueda, and Nic Shannon, “Theory of inelastic neutron scattering in a field-induced spin-nematic state,” Phys. Rev. B 91, 174402 (2015).

[63] S. V. Isakov, R. Moessner, and S. L. Sondhi, “Why spin ice obeys the ice rules,” Phys. Rev. Lett. 95, 217201 (2005).

[64] Owen Benton, L. D. C. Jaubert, Han Yan, and Nic Shannon, “A spin-liquid with pinch-line singularities on the pyrochlore lattice,” Nat Commun 7, 11572 (2016).

[65] P. A. McClarty, S. H. Curnoe, and M. J. P. Gingras, “Energetic selection of ordered states in a model of the Er₂Ti₂O₇ frustrated pyrochlore XY antiferromagnet,” Journal of Physics: Conference Series 145, 012032 (2009).

[66] H. Yan, O. Benton, L. Jaubert, and N. Shannon, “Living on the edge: ground-state selection in quantum spin-ice pyrochlores,” arXiv:1311.3501.

[67] L. D. C. Jaubert, Owen Benton, Jeffrey G. Rau, J. Oitmaa, R. R. P. Singh, Nic Shannon, and Michel J. P. Gingras, “Are Multiphase Competition and Order by Disorder the Keys to Understanding Yb₂Ti₂O₇?” Phys. Rev. Lett. 115, 267208 (2015).

[68] S. E. Palmer and J. T. Chalker, “Order induced by dipolar interactions in a geometrically frustrated antiferromagnet,” Phys. Rev. B 62, 488–492 (2000).

[69] Mathieu Taillefumier, Julien Robert, Christopher L. Henley, Roderich Moessner, and Benjamin Canals, “Semiclassical spin dynamics of the antiferromagnetic heisenberg model on the kagome lattice,” Phys. Rev. B 90, 064419 (2014).

[70] J. Knolle, D. L. Kovrizhin, J. T. Chalker, and R. Moessner, “Dynamics of a two-dimensional quantum spin liquid: Signatures of emergent majorana fermions and fluxes,” Phys. Rev. Lett. 112, 207203 (2014).

[71] Matthias Punk, Debanjan Chowdhury, and Subir Sachdev, “Topological excitations and the dynamic structure factor of spin liquids on the kagome lattice,” Nat Phys 10, 289–293 (2014).

[72] Samuel Bieri, Laura Messio, Bernard Bernu, and Claire Lhuillier, “Gapless chiral spin liquid in a kagome heisenberg model,” Phys. Rev. B 92, 060407 (2015).

[73] Tian-Heng Han, Joel S. Helton, Shaoyan Chu, Daniel G. Nocera, Jose A. Rodriguez-Rivera, Collin Broholm, and Young S. Lee, “Fractionalized excitations in the spin-liquid state of a kagome-lattice antiferromagnet,” Nature 492, 406–410 (2012).

[74] Joseph A. M. Paddison, Marcus Daum, Zhiling Dun, Georg Ehlers, Yaohua Liu, Matthew B. Stone, Haidong Zhou, and Martin Mourigal, “Continuous excitations of the triangular-lattice quantum spin liquid YbMgGaO₄,” Nat Phys 13, 117–122 (2017).

[75] Lucile Savary and Leon Balents, “Quantum spin liquids: a review,” Reports on Progress in Physics 80, 016502 (2017).

[76] Lucile Savary and T. Senthil, “Probing Hidden Orders with Resonant Inelastic X-Ray Scattering,” arXiv:1506.04752.

[77] Andrew Smerald and Nic Shannon, “Theory of spin excitations in a quantum spin-nematic state,” Phys. Rev. B 88, 184430 (2013).

[78] Oleg A. Starykh and Leon Balents, “Excitations and quasi-one-dimensionality in field-induced nematic and spin density wave states,” Phys. Rev. B 89, 104407.
(2014).

[79] Andrew Smerald and Nic Shannon, “Theory of NMR $1/T_1$ relaxation in a quantum spin nematic in an applied magnetic field,” Phys. Rev. B **93**, 184419 (2016).

[80] Solène Guitteny, Julien Robert, Pierre Bonville, Jacques Ollivier, Claudia Decorse, Paul Steffens, Martin Boehm, Hannu Mutka, Isabelle Mirebeau, and Sylvain Petit, “Anisotropic Propagating Excitations and Quadrupolar Effects in $Tb_2Ti_2O_7$,” Phys. Rev. Lett. **111**, 087201 (2013).

[81] Yao-Dong Li and Gang Chen, “Symmetry enriched U(1) topological orders for dipole-octupole doublets on a pyrochlore lattice,” Phys. Rev. B **95**, 041106 (2017).

[82] Z. L. Dun, M. Lee, E. S. Choi, A. M. Hallas, C. R. Wiebe, J. S. Gardner, E. Arrighi, R. S. Freitas, A. M. Arevalo-Lopez, J. P. Attfield, H. D. Zhou, and J. G. Cheng, “Chemical pressure effects on magnetism in the quantum spin liquid candidates $Yb_2X_2O_7$ ($X = Sn, Ti, Ge$),” Phys. Rev. B **89**, 064401 (2014).

[83] C. R. Wiebe and A. M. Hallas, “Frustration under pressure: Exotic magnetism in new pyrochlore oxides,” APL Materials **3**, 041519 (2015).

[84] K Matsuhira, C Sekine, C Paulsen, and Y Hinatsu, “Low-temperature magnetic properties of the geometrically frustrated pyrochlore $Pr_2Sn_2O_7$,” Journal of Magnetism and Magnetic Materials **272–276**, E981 – E982 (2004).

[85] J. Robert, E. Lhotel, G. Remenyi, S. Sahling, I. Mirebeau, C. Decorse, B. Canals, and S. Petit, “Spin dynamics in the presence of competing ferromagnetic and antiferromagnetic correlations in $Yb_2Ti_2O_7$,” Phys. Rev. B **92**, 064425 (2015).

[86] J. D. Thompson, P. A. McClarty, D. Prabhakaran, I. Cabrera, T. Guidi, and R. Coldea, “Quasiparticle Breakdown and Spin Hamiltonian of the Frustrated Quantum Pyrochlore $Yb_2Ti_2O_7$ in Magnetic Field,” arXiv:1703.04506.

[87] K. A. Ross, J. P. C. Ruff, C. P. Adams, J. S. Gardner, H. A. Dabkowska, Y. Qiu, J. R. D. Copley, and B. D. Gaulin, “Two-Dimensional Kagome Correlations and Field Induced Order in the Ferromagnetic XY Pyrochlore $Yb_2Ti_2O_7$,” Phys. Rev. Lett. **103**, 227202 (2009).

[88] A. Maisuradze, P. Dalmas de Réotier, A. Yaouanc, A. Forget, C. Baines, and P. J. C. King, “Anomalously

slow spin dynamics and short-range correlations in the quantum spin ice systems $Yb_2Ti_2O_7$ and $Yb_2Sn_2O_7$,” Phys. Rev. B **92**, 094424 (2015).

[89] Yukio Yasui, Minoru Soda, Satoshi Iikubo, Masafumi Ito, Masatoshi Sato, Nobuko Hamaguchi, Taku Matsushita, Nobuo Wada, Tetsuya Takeuchi, Naofumi Aso, and Kazuhisa Kakurai, “Ferromagnetic Transition of Pyrochlore Compound $Yb_2Ti_2O_7$,” Journal of the Physical Society of Japan **72**, 3014–3015 (2003).

[90] E. Lhotel, S. R. Giblin, M. R. Lees, G. Balakrishnan, L. J. Chang, and Y. Yasui, “First-order magnetic transition in $Yb_2Ti_2O_7$,” Phys. Rev. B **89**, 224419 (2014).

[91] A. M. Hallas, J. Gaudet, M. N. Wilson, T. J. Munsie, A. A. Aczel, M. B. Stone, R. S. Freitas, A. M. Arevalo-Lopez, J. P. Attfield, M. Tachibana, C. R. Wiebe, G. M. Luke, and B. D. Gaulin, “XY antiferromagnetic ground state in the effective $S = \frac{1}{2}$ pyrochlore $Yb_2Ge_2O_7$,” Phys. Rev. B **93**, 104405 (2016).

[92] Stephen Powell, “Ferromagnetic coulomb phase in classical spin ice,” Phys. Rev. B **91**, 094431 (2015).

[93] M. E. Brooks-Bartlett, S. T. Banks, L. D C Jaubert, A. Harman-Clarke, and P. C W Holdsworth, “Magnetic-Moment Fragmentation and Monopole Crystallization,” Phys. Rev. X **4**, 011007 (2014).

[94] Owen Benton, “Quantum origins of moment fragmentation in $Nd_2Zr_2O_7$,” Phys. Rev. B **94**, 104430 (2016).

[95] E. Lefrançois, V. Cathelin, E. Lhotel, J. Robert, P. Lejay, C. V. Colin, B. Canals, F. Damay, J. Ollivier, B. Fåk, L. C. Chapon, R. Ballou, and V. Simonet, “Magnetic charge injection in spin ice: a new way to fragmentation,” arXiv:1702.02864.

[96] K. A. Ross, J. W. Krizan, J. A. Rodriguez-Rivera, R. J. Cava, and C. L. Broholm, “Static and dynamic $xy$-like short-range order in a frustrated magnet with exchange disorder,” Phys. Rev. B **93**, 014433 (2016).

[97] K. A. Ross, J. M. Brown, R. J. Cava, J. W. Krizan, S. E. Nagler, J. A. Rodriguez-Rivera, and M. B. Stone, “Single-ion properties of the $S_{eff} = \frac{1}{2}$ XY antiferromagnetic pyrochlores $NaA'Co_2F_7$ ($A' = Ca^{2+}, Sr^{2+}$),” Phys. Rev. B **95**, 144414 (2017).

[98] Hirokazu Tsunetsugu, “Antiferromagnetic Quantum Spins on the Pyrochlore Lattice,” Journal of the Physical Society of Japan **70**, 640–643 (2001).