Quantum versus Classical Spin Fragmentation in Dipolar Kagome Ice Ho3Mg2Sb3O14

A promising route to realize entangled magnetic states combines geometrical frustration with quantum-tunneling effects. Spin-ice materials are canonical examples of frustration, and Ising spins in a transverse magnetic field are the simplest many-body model of quantum tunneling. Here, we show that the tripod kagome lattice material Ho${_3}$Mg${_2}$Sb${_3}$O${_{14}}$ unites an ice-like magnetic degeneracy with quantum-tunneling terms generated by an intrinsic splitting of the Ho$^{3+}$ ground-state doublet, which is further coupled to a nuclear spin bath. Using neutron scattering and thermodynamic experiments, we observe a symmetry-breaking transition at 0.32 K to a remarkable quantum state with three peculiarities: a dramatic recovery of magnetic entropy associated with the strongly coupled electronic and nuclear spins; a fragmentation of the spin into periodic and ice-like components strongly affected by quantum fluctuations; and persistent inelastic magnetic excitations down to 0.12 K. These observations deviate from expectations of classical spin fragmentation physics on a kagome lattice, but can be understood within a framework of dipolar kagome ice under a homogeneous transverse magnetic field. By combining our experimental data with exact diagonalization and mean-field modeling, we establish the existence of a highly entangled fragmented state in a region where the transverse field remains a perturbation to the dipolar interactions, which we coin as a quantum spin fragmented. However, hyperfine interactions play a crucial role in suppressing quantum correlations and dramatically alter the single-ion and collective properties of Ho${_3}$Mg${_2}$Sb${_3}$O${_{14}}$. Our results highlight the crucial role played by hyperfine interactions in frustrated quantum magnets, and motivate further theoretical investigations of spin fragmentation and coherent quantum tunneling.

A promising route to realize entangled magnetic states combines geometrical frustration with quantumtunneling effects. Spin-ice materials are canonical examples of frustration, and Ising spins in a transverse magnetic field are the simplest many-body model of quantum tunneling. Here, we show that the tripod kagome lattice material Ho3Mg2Sb3O14 unites an ice-like magnetic degeneracy with quantum-tunneling terms generated by an intrinsic splitting of the Ho 3+ ground-state doublet, which is further coupled to a nuclear spin bath. Using neutron scattering and thermodynamic experiments, we observe a symmetry-breaking transition at T * ≈ 0.32 K to a remarkable quantum state with three peculiarities: a dramatic recovery of magnetic entropy associated with the strongly coupled electronic and nuclear degrees of freedom; a fragmentation of the spin into periodic and ice-like components strongly affected by quantum fluctuations; and persistent inelastic magnetic excitations down to T ≈ 0.12 K. These observations deviate from expectations of classical spin fragmentation physics on a kagome lattice, but can be understood within a framework of dipolar kagome ice under a homogeneous transverse magnetic field. By combining our experimental data with exact diagonalization and mean-field modeling, we establish the existence of a highly entangled fragmented state in a region where the transverse field remains a perturbation to the dipole-dipole interactions, which we coin as a quantum spin fragmented state. However, hyperfine interactions play a crucial role in suppressing quantum correlations and dramatically alter the single-ion and collective properties of Ho3Mg2Sb3O14. Our results highlight the crucial role played by hyperfine interactions in frustrated quantum magnets, and motivate further theoretical investigations of the interplay between spin fragmentation and coherent quantum tunneling.

I. INTRODUCTION
Quantum spin liquids are exotic states of magnetic matter in which conventional magnetic order is suppressed by strong quantum fluctuations [1]. Frustrated magnetic materials, which have a large degeneracy of classical magnetic ground states, are often good candidates to search for this elusive behavior. A canonical example of frustration is spin ice, in which Ising spins occupy a pyrochlore lattice of corner-sharing tetrahedra [2,3]. Classical ground states obey the "two in, two out" ice rule for spins on each tetrahedron, and thermal excitations behave as deconfined magnetic monopoles [4][5][6][7]. These pairs of fractionalized excitations interact via Coulomb's law and correspond to topological defects of a classical field theory obtained by coarsegraining spins into a continuous magnetization. In principle, topological quantum excitations can be generated by adding quantum-tunneling terms to the classical spin-ice model-Z. Dun and X. Bai contributed equally to this work. * zdun3@gatech.edu † paddisonja@ornl.gov ‡ mourigal@gatech.edu § hzhou10@utk.edu e.g., by adding couplings between the transverse components of spins [8][9][10], or by introducing a local magnetic field transverse to the Ising spins [11][12][13][14]. A search for real materials that realize such quantum spin-ice states has found several promising candidates (see, e.g., [15][16][17][18][19][20][21][22][23][24][25]). However, important challenges remain, including the determination of the often-complex spin Hamiltonian [26][27][28], the subtle role that structural disorder may play [29][30][31], and the computational challenges associated with simulations of three-dimensional (3D) quantum magnets [32,33].
A promising alternative route towards quantum analogs of spin ice is offered by two-dimensional Ising ferromagnets on a kagome lattice. When the spins are confined to point either towards or away from the center of each triangle of the lattice, a highly degenerate kagome ice state is stabilized with a "one in, two out" or "two in, one out" local ice rule on each triangle [35]. Both a quantum-tunneling term and an external magnetic field are required to enable tunneling between these states [36,37]. Remarkably, when the long range magnetic dipole-dipole interaction is introduced, the effective Coulomb interaction between emergent magnetic charges-defined in Fig. 1(a)-selects a sub-space of the kagome ice manifold and drives a phase transition at low temperatures to a state with staggered emergent-charge ordering [38,39]. Nevertheless, this state still possesses nonzero entropy, because each emer-  34. The expectation values of spins σ z i are represented by black arrows. Each triangle has one spin pointing "in" (towards its center) and two pointing "out" (away from its center), or vice versa. The emergent magnetic charge of a triangle is defined as the number of spins pointing "in" minus the number pointing "out". Positive (Qj = +) and negative (Qj = −) emergent charges are represented as red and blue triangles, respectively, and form a staggered arrangement. Three distinct spin configurations are possible for a given emergent charge, which yields a macroscopic number of degenerate spin configurations associate with emergent charge ordering. Spin fragmentation decomposes each unit-length spin into "divergence-full" and "divergence-free" channels (center and right images, respectively). The fragmented spins are shown as orange circles with diameter proportional to the length of the fragmented spin. The green hexagon represents the flipping of six spins around a closed loop: this is the simplest process that connects two distinct spin configurations within the degenerate CSF manifold. (b) Conceptual zero temperature phase diagram of a dipolar kagome ice model under a transverse field, as informed by our exact diagonalization on small clusters (Section VI). The strength of the transverse field h x relative to the scale of the dipolar interaction D yields three distinct regimes: (i) h x D (perturbative regime) where quantum dynamics is generated by tunneling processes on closed hexagonal loops, a regime we call "quantum spin fragmented" (QSF); (ii) h x 0.2D (non-perturbative regime) where collective quantum dynamics persist but the transverse field leads to the disappearance of emergent charges by mixing states from within and outside the CSF manifold; (iii) h x D, a crossover to the single-ion regime.
gent charge retains a threefold degeneracy of spin orientations [38]; hence, ordering of the emergent charges does not imply complete ordering of the spins. Fig. 1(a) shows that such spin structures can be decomposed into a "divergence-full" channel in which spins are spatially ordered, and an "divergence-free" channel in which spins remain spatially disordered-a process known as spin fragmentation [34,40]. Neutron-scattering measurements provide a direct experimental signature of spin fragmentation via the coexistence of magnetic Bragg peaks and highly-structured magnetic diffuse scattering with pinchpoint singularities [21,41,42]. The divergence-full channel corresponds to an "all-in/all-out" (AIAO) order of fragmented spins that reflects the long-range staggered arrangement of emergent charges. In the divergence-free channel, fragmented spins are disordered but correlated and the constraint that every triangle has zero emergent charge yields a Coulomb phase analogous to pyrochlore spin ices [40]. Without additional quantum effects, the classical spinfragmented (CSF) state described above can be viewed as a classical spin-liquid coexisting with magnetic order. Its elementary excitations are thermally flipped spins that map onto pairs of magnetic monopoles, i.e., defects in the divergencefree channel. Typically, such thermally-activated excitations are exponentially suppressed at low temperature, as observed experimentally in the CSF kagome ice compound Dy 3 Mg 2 Sb 3 O 14 [41,43]. The question arises whether it is possible to stabilize a quantum state with collective dynamics revealing fluctuations between degenerate classical spin configurations, in close analogy to the quantum ice physics proposed for pyrochlore and square-lattice systems [8-10, 12, 14]. Conceptually, a magnetic field transverse to the Ising moments can yield quantum-tunneling processes at low temperature that connect different CSF states through the concurrent flipping of six spins, a process shown with a green hexagon in Fig. 1(a). In this scenario, we still expect the coexistence of Bragg peaks and highly-structured diffuse scattering in magnetic neutron scattering experiments. However, because magnetic correlations are no longer static, we expect this state to host coherent collective excitations at low temperatures, akin to the emergent monopole and photon-like excitations in quantum spin ice [10,44]. To distinguish it from the CSF, we tentatively describe this putative new state as "quantum spin fragmented" (QSF).
In this work, we present comprehensive inelastic neutron scattering data which shows that spin dynamics exist at the lowest measurable temperatures (≈ 0.1 K) in the dipolar kagome Ising magnet Ho 3 Mg 2 Sb 3 O 14 [45]. This material is one of a series of "tripod kagome" materials derived from the pyrochlore structure by chemical substitution, yielding kagome planes of magnetic rare-earth ions separated by triangular planes of nonmagnetic Mg 2+ ions [ Fig. 2(a)]. Previous measurements of isostructural Dy 3 Mg 2 Sb 3 O 14 revealed a CSF state at low temperature [41], in which no spin dynamics were observed in either neutron-scattering or ac-susceptibility data [41,46]. Our measurements on Ho 3 Mg 2 Sb 3 O 14 uncover a spin-fragmented state with instantaneous magnetic correlations closely resembling that of Dy 3 Mg 2 Sb 3 O 14 . Yet, we observe structured dynamic magnetic correlations at low temperature, indicating persistent spin dynamics in sharp contrast to the Dy 3+ compound. We show this stems from the low symmetry of the tripod-kagome structure and the non-Kramers nature of the Ho 3+ ion, the combination of which generates an effective local magnetic field transverse to the Ising magnetic dipole moments. The effective low-energy Hamiltonian for Ho 3 Mg 2 Sb 3 O 14 thus maps onto an iconic model of quantum magnetism-interacting Ising spins in a transverse magnetic field [47]. We use neutron-scattering experiments to determine this Hamiltonian, and employ a combination of exact diagonalization, field theoretic, and Monte Carlo methods to understand its spin correlations and excitations. Remarkably, our calculations suggest that a QSF state is realized in models of dipolar kagome ice in a small transverse field. In Ho 3 Mg 2 Sb 3 O 14 , however, the single-ion physics of Ho 3+ ions is profoundly affected by the strong nuclear hyperfine coupling, which eventually destroys most coherent quantum effects between sites.
Our paper is structured as follows.
In Section II, we summarize the experimental methods that we employ. In Section III, we present neutron-scattering measurements of the crystal-field excitations of Ho 3 Mg 2 Sb 3 O 14 and of a structurally analogous but magnetically-dilute system (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 , and use these measurements to parameterize the spin Hamiltonian of Ho 3 Mg 2 Sb 3 O 14 . In Section IV, we report heat-capacity measurements that identify a magnetic phase transition at T * = 0.32 K accompanied by a large specific heat feature of coupled electronic-nuclear origin. In Section V, we report low-temperature inelastic neutron-scattering measurements on polycrystalline samples of Ho 3 Mg 2 Sb 3 O 14 . They reveal that spin fragmentation occurs below T * , and that low-energy spin excitations are structured in both momentum and energy space at the lowest measurable temperatures. In Section VI, we use theoretical modeling to understand our data. Finally, we conclude in Section VII with a discussion of the general implications of our study.

II. METHODS
Two different polycrystalline samples of Ho 3 Mg 2 Sb 3 O 14 were prepared for this study, the first using a traditional solid state reaction method (referred to as s.s. sample), and the second using a sol-gel method (referred to as s.g. sample). For the s.s. sample, stoichiometric ratios of Ho 2 O 3 (99.9%), MgO (99.99%), and Sb 2 O 3 (99.99%) fine powder were carefully ground and reacted at a temperature of 1350 • C in air for 24 hours. This heating step was repeated until the amount of impurity phases as determined by X-ray diffraction was not reduced further. The synthesized s.s. sample contained a small amount of Ho 3 SbO 7 impurity (2.29(18) wt%), which orders antiferromagnetically at T N = 2.07 K [48]. This impurity can be removed by the sol-gel synthesis method. For this synthesis, stoichiometric amounts of Ho(NO 3 ) 3 , Mg(NO 3 ) 3 (prepared by dissolving Ho 2 O 3 and MgO in hot diluted nitric acid solution), and antimony tartarate (prepared by dissolving Sb 2 O 3 in hot tartaric acid solution) were first mixed in a beaker. Citrate acid with a metal-to-citrate molar ratio of 1:2 was then added to the solution followed by a subsequent heating on a hot plate at 120 • C overnight to remove excessive water. The obtained gel-like solution was slowly heated to 200 • C in a box furnace to decompose the nitrate, and was pyrolyzed at 600 • C for 10-12 h in air. The obtained powder was then ground, pressed into a pellet and re-heated at 1300 • C until a well-reacted crystalline powder was obtained. Heat-capacity measurements presented below show that the thermo-magnetic behavior of the two samples is almost identical, expect for a small peak around 2.1 K in the s.s. sample originating from the Ho 3 SbO 7 impurity. A Ho-diluted (Ladoped) sample of (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 was also synthesized with the same sol-gel technique, with 99% Ho 2 O 3 replaced by La 2 O 3 powder (99.9%, baked at 900 • C overnight before use).
Low-temperature specific-heat measurements were performed on a Quantum Design Physical Properties Measurement System instrument using dilution refrigerator (0.07 ≤ T ≤ 4K) and standard (1.6 ≤ T ≤ 100K) probes. For the dilution refrigerator measurement, the powder samples were cold-sintered with Ag powder. The contribution of the Ag powder was measured separately and subtracted from the data. The lattice contribution to the heat capacity was estimated from measurements of the isostructural nonmagnetic compound La 3 Mg 2 Sb 3 O 14 .
Powder X-ray diffraction measurements were carried out with Cu Kα radiation (λ = 1.5418Å) in transmission mode. Powder neutron-diffraction measurements were carried out using the HB-2A high-resolution powder diffractometer [49] at the High Flux Isotope Reactor at Oak Ridge National Laboratory, with a neutron wavelength of 1.546Å. Rietveld refinements of the crystal and magnetic structures were carried out using the FULLPROF suite of programs [50]. Peakshapes were modeled by Thompson-Cox-Hastings pseudo-Voigt functions, and backgrounds were fitted using Chebyshev polynomial functions.
Inelastic neutron-scattering measurements on the s.s. sample of Ho 3 Mg 2 Sb 3 O 14 were carried out using the Fine-Resolution Fermi Chopper Spectrometer (SEQUOIA) [51] at the Spallation Neutron Source of Oak Ridge National Laboratory, and the Disk Chopper Spectrometer (DCS) [52] at the NIST Center for Neutron Research. For the SEQUOIA experiment, a ∼5 g powder sample of Ho 3 Mg 2 Sb 3 O 14 was loaded in an aluminum sample container and cooled to 4 K with a closed-cycle refrigerator. Data were measured with incident neutron energies of 120, 60, and 8 meV. The same measurements were repeated for an empty aluminum sample holder and used for background subtraction. For the DCS measurements, the same sample was loaded in a copper can, filled with 10 bar of helium gas at room temperature, and cooled to millikelvin temperatures using a dilution refrigerator. The measurements were carried out with an incident neutron energy of 3.27 meV at temperatures between 0.12 and 40 K. Measurements of an empty copper sample holder were also made and used for background subtractions. Due to the large specific heat and related relaxation processes below 1 K, a thermal stabilization time of 6 h was used; no change in the data was observed after this waiting time. Data reduction was performed using the DAVE program [53]. For modeling and fitting purposes, data were corrected for background scattering using empty-container measurements and/or high-temperature measurements, as specified in the text. These data were also corrected for neutron absorption [54], placed on an absolute intensity scale by scaling to the nuclear Bragg profile, and the magnetic scattering from the Ho 3 SbO 7 impurity below its T N of 2.07 K was subtracted as described in Ref. 41. Additional inelastic neutron-scattering measurements on Ho 3 Mg 2 Sb 3 O 14 (s.g. sample) and (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 were carried out using the OSIRIS backscattering spectrometer at the ISIS neutron source with a final neutron energy of 1.84 meV. For (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 , a 13.6 g powder sample was loaded into an aluminum can and was cooled with an Orange cryostat to the base temperature of 1.6 K. For Ho 3 Mg 2 Sb 3 O 14 , 4 g of s.g. powder was wrapped in a thin copper foil and placed into a copper sample can with a capillary that allowed helium filling at low temperature. The system was cooled using a dilution refrigerator, and with a maximum of 3 bar filled helium, the lowest sample temperature accessed in this experiment was estimated to be 400±50 mK by comparison with the DCS data.
For convenience, in the following sections, we use a unit system with k B = 1 and = 1, so that all energies are given in units of K.

A. Crystal structure and interactions
The crystal structure of Ho 3 Mg 2 Sb 3 O 14 (space group R3m) is shown in Fig. 2(a), and contains kagome planes of magnetic Ho 3+ ions separated by triangular layers of nonmagnetic Mg 2+ [45]. The Ho 3+ site has C 2h point symmetry and its local environment contains eight oxygen atoms [45,46]. The orientations of Ho 3+ magnetic dipole moments are constrained by crystal electric field effects to point along the line connecting Ho 3+ to its two closest oxygen neighbors, which are situated near the centroids of the MgHo 3 tetrahedra [ Fig. 2(b)]. Rietveld co-refinements to X-ray and neutron powder-diffraction data for Ho 3 Mg 2 Sb 3 O 14 (s.s. sample) confirm this crystal structure, and reveal a small amount of Ho 3+ /Mg 2+ site mixing such that 3.2(2)% of Ho 3+ atomic positions are occupied in a disordered way by Mg 2+ (see Appendix A). Hence, the extent of chemical disorder in Ho 3 Mg 2 Sb 3 O 14 is less than in its Dy 3+ analog, where the corresponding value is 6(2)% [41].
We anticipate that the spin Hamiltonian for Ho 3 Mg 2 Sb 3 O 14 may be written as a sum of three terms, where H cf , H hf , and H int denote respectively the crystalfield, nuclear hyperfine, and pairwise interaction Hamiltonian. We now consider the origin, form, and magnitude of each term, and show that they are all relevant in Ho 3 Mg 2 Sb 3 O 14 at low temperature.

B. Crystal-field Hamiltonian
We use inelastic neutron-scattering measurements and point-charge calculations to determine the parameters of H cf . The high-energy spectrum observed in Ho 3 Mg 2 Sb 3 O 14 (s.s. sample) comprises five crystal-field excitations, with energies and relative intensities that resemble those of pyrochlore spin ice Ho 2 Ti 2 O 7 [55,56] except for an overall downwards renormalization in energy [ Fig. 2(c)]. This overall resemblance is expected given the similar local environments for Ho 3+ ions in these two systems. However, a crucial difference stems from the reduced C 2h symmetry of the Ho 3+ site in Ho 3 Mg 2 Sb 3 O 14 compared to the D 3d symmetry in Ho 2 Ti 2 O 7 . Whereas the crystal-field ground-state in Ho 2 Ti 2 O 7 is a non-Kramers doublet, in Ho 3 Mg 2 Sb 3 O 14 all crystal-field levels are necessarily singlets [45]. However, as a probable consequence of spin-spin interactions (see Section V), the excitation associated with the splitting of the groundstate doublet is strongly overdamped [ Fig. 2(f)]. This precludes a direct neutron-scattering measurement of this energy splitting in Ho 3 Mg 2 Sb 3 O 14 . We addressed this problem using two complementary approaches. First, we performed point-charge calculations using an effective charge model that matches the high-energy crystal-field excitation spectrum (see Appendix B). Second, we validated the low-energy predictions of this model using high-resolution neutron scattering measurements on a magnetically-dilute Ho tripod kagome compound, (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 , such that interaction effects between sites are negligible. We discuss these results in turn below.
Our point-charge model predicts that the ground-state doublet is split into two singlets separated by an energy gap ∆ ≈ 1.74 K (see Appendix B). The two singlets are well approximated by symmetric and anti-symmetric superposition of pure free-ion states, where |± ≡ |J = 8, J z = ±8 represents the ground-state doublet of Ho 3+ in Ho 2 Ti 2 O 7 if we ignore the small contributions from other J z components [55,57]. At low temperatures, only |0 and |1 are thermally populated because of their 190 K separation from higher-energy crystal-field levels [ Fig. 2(c)]. The form of |0 and |1 allows for a nonzero angular momentum matrix element 0|Ĵ z |1 = 7.79 while 0|Ĵ α |1 with α = x, y are vanishingly small (see Appendix B). Similar to the procedure for pyrochlore spin ices [57], we construct Pauli matrices using the two thermally-accessible crystal-field states, In this framework, the total magnetic dipole moment operator is related to σ z as where g J = 5 4 is the Ho 3+ Landé factor, andẑ is a local Ising axis shown in Fig. 2(b). The nonzero matrix element 0|Ĵ z |1 is therefore expected to generate a total magnetic moment of magnitude 9.74 µ B , which is conserved. In contrast, static moments only appear when σ z becomes non-zero, e.g. under an external magnetic field. This prediction is supported by our isothermal magnetization measurements between 1.8 K and 40 K, which are consistent with a model of paramagnetic Ising spins (see Appendix C).
It is established [14,47] that an energy splitting between two crystal-field singlets can be exactly mapped into a transverse magnetic field acting on a corresponding doublet. This  (2), 238(2), 388(9), 627(5), and 743(9) K, whose peak positions are extracted from Lorentzian fits to the data (black dashed lines). The overall fit, including a modeled phonon background, is shown as a red line. (d) Comparison of the single-ion response expected for a Ho 3+ ion in a pyrochlore vs. a tripod kagome geometry. Crystal-field singlets |0 and |1 in Ho3Mg2Sb3O14 can be effectively described as a transverse field h x acting on the non-Kramers doublet |± of Ho2Ti2O7. While hyperfine interactions split the electronic-nuclear manifold into eight levels in both systems, the transverse field h x transforms the classical single-ion picture of Ho2Ti2O7 into a quantum picture in Ho3Mg2Sb3O14, where each eigenstate is a quantum superposition of |+, I z and |−, I z . Red dashed arrows mark the allowed transitions between eigenstates of the same I z . (e) Calculated scattering intensity of single-ion excitations as a function of h x at a temperature of 1.6 K. Dashed lines present the transition energies whose intensities are calculated and convoluted with a Lorentzian function with FWHM = 0.35 K. (f) Low energy inelastic neutron scattering response for pure Ho3Mg2Sb3O14 (s.g. sample) and the dilute Ho-tripod magnet (Ho0.01La0.99)3Mg2Sb3O14 at a temperature of 1.6 K measured at OSIRIS. The two datasets are normalized to the same total spectrum weight. Solid lines represent the single-ion spectrum calculated using h x = 0.57 K, which yields the best agreement with the experimental data. mapping can be understood by recognizing that |0 and |1 are the eigenstates of the σ x = 0 1 1 0 Pauli matrix. Therefore, in the Pauli matrix picture, our crystal field Hamiltonian can be recast as where h x = ∆/2 is the intrinsic transverse field [ Fig. 2(d)].
We will use this pseudo-spin representation in the rest of this paper.

C. Nuclear hyperfine Hamiltonian
Hyperfine interactions couple non-zero nuclear spins to the local magnetic field from surrounding electrons. With a nuclear spin quantum number I = 7/2, 165 Ho is the only stable isotope of holmium and its hyperfine energy-scale is the strongest among the rare-earth elements. For a non-Kramers electronic system in the pseudo-spin approximation, the hyperfine Hamiltonian takes the simple form [58] where I z = −I, ..., I labels the z-component of the nuclear spin operator, and A hf = 0.319 K is the hyperfine coupling constant for Ho [59,60]. We neglect the electric quadrupole coupling constant because its energy scale (P = 0.004 K) is very small. In Ho 2 Ti 2 O 7 , the hyperfine coupling splits the combined electronic and nuclear system into eight uniformly-spaced levels; each level is doubly degenerate because of the Kramers degeneracy of the combined electronic and nuclear spin system [ Fig. 2(d)]. At the single-ion level, the system remains classical because both I z and σ z are good quantum numbers, and eigenstates can be labeled as |±, I z . By contrast, for Ho 3 Mg 2 Sb 3 O 14 , the single-ion Hamiltonian contains both an intrinsic transverse field and hyperfine interactions, This generates a dynamic electronic moment since neither σ z nor σ x are good quantum numbers. Diagonalizing the 16-dimensional H SI yields an unevenly spaced spectrum for which each eigenstate is a quantum superposition of |+, I z and |−, I z whose mixing depends on the values of I z and h x , and the lowest energy states always have I z = ±7/2. To illustrate this effect, Fig. 2(e) shows the low-temperature magnetic scattering intensity of electronic spins as a function of h x . For vanishing h x , the magnetic response is purely elastic because the dipolar matrix elements connecting different eigenstates are zero. With increasing h x , the scattering acquires an inelastic component that comprises four excitations, arising from transitions between eigenstates with the same I z [ Fig. 2 To test our model single-ion Hamiltonian against experiment, we employ high-resolution neutron scattering measurements of the magnetically-dilute Ho tripod kagome compound (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 . We assume that the Ho 3+ ions are randomly distributed on the kagome lattice and pairwise interactions between them are negligible. Neutron-scattering data measured at 1.6 K are shown in Fig. 2(f); they display a broad and intense peak at ω = 1.3 K, along with two narrower and weaker side peaks at ω = 1.9 K and 2.5 K. The experimental data are in excellent quantitative agreement with exact diagonalization calculations of H SI , taking . This result yields strong evidence that H SI describes well the single-ion properties of Ho-based tripod kagome magnets. Importantly, due to the larger ionic size of La 3+ compared with Ho 3+ , the lattice parameters of (Ho 0.01 La 0.99 ) 3 Mg 2 Sb 3 O 14 are approximately 3% larger than those of Ho 3 Mg 2 Sb 3 O 14 . Using a power law scaling of h x as a function of lattice parameters, we extrapolate to h x = 0.85 K in Ho 3 Mg 2 Sb 3 O 14 , which is in good agreement with the point charge estimate of h x = 0.77 K. We take h x = 0.85 K for Ho 3 Mg 2 Sb 3 O 14 throughout the rest of this paper.

D. Transverse Ising model
We now consider the effect of pairwise interactions J ij between Ho 3+ ions. By analogy with spin-ice pyrochlores [3] and isostructural Dy 3 Mg 2 Sb 3 O 14 [41], we expect that J ij contains a combination of nearest-neighbor exchange interactions J nn and long-range magnetic dipolar interactions of overall scale D. In principle, interactions between transverse spin components are also possible, but they are expected to be several orders of magnitude smaller [57], so we do not consider them further. The pairwise interaction Hamiltonian is therefore with where δ rij ,rnn is the Kronecker delta function, r nn is the distance between nearest-neighbor Ho 3+ ions, r ij is the distance between ions at positions r i and r j , andr ij = (r i − r j )/r ij . The value of D = µ 0 µ 2 /(4πk B r 3 nn ) = 1.29 K is fixed by the crystal structure, and we will obtain an experimental estimate of J nn ≈ −0.64 K in Section V. Withẑ i ·ẑ j = −0.28, and −3(ẑ i ·r nn )(ẑ j ·r nn ) = 1.93, the sum of exchange and dipolar couplings at the nearest neighbor level is approximately 1.69 K, and hence antiferromagnetic in the pseudo-spin language. Consequently, magnetic interactions between sites are frustrated, and furthermore, comparable in magnitude to h x .
Using the results of the previous subsections, we rewrite the full spin Hamiltonian, Eq. (1), as Eq. (10) is equivalent in form to an Ising model in a transverse field (TIM), with an additional on-site longitudinal field due to the hyperfine coupling. The TIM has been used to model diverse physical phenomena, including ferroelectricity [61,62], superconductivity [63], quantum information [64,65], and quantum phase transitions [66,67]. Typically, the pairwise interactions that drive magnetic ordering compete with the transverse field that drives quantum tunneling. In the absence of geometrical frustration, a phase transition only occurs to a magnetically-ordered state if J ij dominates over h x . The interplay of frustration and transverse field may generate exotic quantum phases [11,14,68,69]. On the kagome lattice, the TIM with nearest-neighbor antiferromagnetic interactions is predicted to have a quantum-disordered ground state for small h x , smoothly connected to a quantum paramagnetic state at large h x [11,68,69]. On the pyrochlore lattice, an external field cannot be applied transverse to all spins simultaneously because the different local Ising axes are not coplanar, and a homogeneous transverse field that emerges at single ion level is absent in chemically-ordered pyrochlores. However, transverse fields generated by spin-spin interactions are related to monopole hoping in spin ice [70]; and transverse fields generated by chemical disorder have been identified as a possible route to pyrochlore QSL states [14], and used to explain the spin dynamics of Pr 2 Zr 2 O 7 [22,24] and Tb 2 Ti 2 O 7 [71,72]. Nevertheless, a potential challenge to modeling such materials is that chemical disorder generates a broad distribution of transverse fields in the sample [73]. Hence, a key feature of Ho 3 Mg 2 Sb 3 O 14 is that its transverse field is intrinsic to the chemically-ordered structure, and to a first approximation is homogeneous.

IV. SPECIFIC-HEAT MEASUREMENTS
We use heat capacity measurements to understand thermodynamic properties of Ho 3 Mg 2 Sb 3 O 14 and identify possible phase transitions. A sharp peak in the magnetic specific heat (C m ) is observed at T * = 0.32 K for both s.s. and s.g. samples, indicating a symmetry-breaking magnetic phase transition [ Fig. 3(a)]. The value of T * is consistent with the broad peak previously observed around 0.4 K using the ac susceptibility technique [45]. Whereas the ac susceptibility peak is frequency dependent [45], the sharpness of the C m peak is inconsistent with a conventional spin freezing scenario. The value of T * is also close to the temperature at which the isostructural compound Dy 3 Mg 2 Sb 3 O 14 undergoes a phase transition from a kagome spin-ice state to a CSF state (∼0.3 K in Ref. [41] and ∼0.37 K in Ref. [46]). As we will show in Section V, T * corresponds to the onset of a spin fragmented state in Ho 3 Mg 2 Sb 3 O 14 , characterized by a reduced ordered moment compared to a CSF state.
Below 1 K, a broad specific heat feature is observed in addition to the sharp peak, consistent with a nuclear spin contribution. By integrating C m /T from 20 K to the lowest measurement temperature of 76 mK, the recovered magnetic entropy reaches 21.5(5) J K −1 mol −1 Ho , which is only 6.5% smaller than the expectation of R(ln 2 + ln 8) = 23.05 J K −1 mol −1 Ho considering both the electronic and nuclear spin degrees of freedom [ Fig. 3(b)]. In Ho-based systems with a doublet singleion ground state, such as Ho metal [59], Ho 2 Ti 2 O 7 [74], and LiHoF 4 [75], the nuclear specific heat consists of a broad peak (Schottky anomaly) below the ordering temperature of the electronic spins. This implies that the dynamics of the electronic and nuclear subsystems separate in such systems, with electronic spins already in their σ z eigenstate when nuclear spins start to follow them at low temperature. This paradigm is not applicable to Ho 3 Mg 2 Sb 3 O 14 because h x mixes |+, I z and |−, I z at each site [ Fig. 2(d)]. Accordingly, the singleion Hamiltonian [Eq. (7)] predicts two peaks for the nuclear specific heat; however, this model strongly disagrees with our experimental data [ Fig. 3(a)]. The observed C m is also very different from that of unfrustrated TIM magnets such as HoF 3 , in which hyperfine interactions act as an effective mean field that precipitates long-range ordering of electronic spins, and the nuclear Schottky anomaly is only manifest below T * [76]. By contrast, in Ho 3 Mg 2 Sb 3 O 14 , more than half of the total magnetic entropy has already been recovered above T * , suggesting the development of short-range spin correlations with a coupled electronic-nuclear character. Our heat capacity measurements thus provide the first experimental hint of many-body physics in Ho 3 Mg 2 Sb 3 O 14 , demonstrating that the hyperfine term, transverse field, and spin-spin interactions in Eq. (10) must be treated on an equal footing.

V. INELASTIC NEUTRON-SCATTERING MEASUREMENTS
We use inelastic neutron-scattering measurements to probe the spin correlations and low-energy spin dynamics of Ho 3 Mg 2 Sb 3 O 14 . Neutron-scattering data as a function of momentum (Q) and energy transfer (ω) are shown in Fig. 4(a) over the temperature range from 4.2 K to 0.12 K. The Qdependence of the magnetic scattering shown in Fig. 4(b) was obtained by integrating over −20 ≤ ω ≤ 20 K and subtracting the paramagnetic data at T = 30 K. Such energy-integrated data measure the Fourier transform of the instantaneous spinpair correlation function. The energy dependence shown in Fig. 4(c) was obtained by integrating the inelastic scattering over 0.4 ≤ Q ≤ 1.6Å −1 and correcting for background scattering using empty-container measurements.
We first discuss the paramagnetic regime above T * . For sample temperatures between 4.2 K and 0.4 K, the energyintegrated magnetic diffuse scattering displays a clear Qdependence, with a broad peak centered at approximately 0.65Å −1 that develops on cooling [ Fig. 4(b)]. This feature closely resembles observations for Dy 3 Mg 2 Sb 3 O 14 , where it was interpreted in terms of the development of kagome-ice correlations with a "one in, two out" or "two in, one out" ice rule on each triangle [41]. The energy-resolved response shows two main magnetic features: an intense and resolutionlimited quasielastic peak; and a broad inelastic tail extending to ω ≈ 15 K that decays slowly with increasing energy transfer [ Fig. 4(c)]. This energy-resolved response differs dramatically from the case of Dy 3 Mg 2 Sb 3 O 14 , in which only elastic (ω = 0) neutron scattering was observed at comparable temperatures, indicating time-independent spin correlations [41]. Moreover, the energy-resolved response of Ho 3 Mg 2 Sb 3 O 14 is very different from the single-ion model discussed in Section III, in which the majority of scattering is expected to be inelastic and concentrated around energy transfer ∆ = 2h x ≈ 1.7 K [ Fig. 2 To obtain a better understanding of our paramagnetic neutron-scattering data, we employ a reciprocal-space meanfield approximation to model the wave-vector dependence of the magnetic diffuse scattering. This approach is exact in the high-temperature limit, and introduces the effect of local spin correlations via a reaction-field term that is determined selfconsistently [83]. The neutron-scattering intensity is calculated via the dynamical susceptibility, which is approximated as where χ 0 (ω) is the single-ion susceptibility obtained by exact diagonalization of H SI , λ is the reaction field [83,84], and µ ∈ {1, 3} labels the normal modes of the Ising system in the same way as for classical mean-field theories [85]. Full details of this method are given in Ref. 83 and references therein. We employ this approach to fit the value of J nn to the energyintegrated magnetic diffuse scattering. The value of J nn is the only free parameter because h x = 0.85 K is fixed (see Section III B). Our fits yield good agreement with the energyintegrated experimental data at 4.2 K [ Fig. 4(b)] and at 10 K (not shown). However, because damping effects are not included in the mean-field calculation, the energy-resolved response deviates from the experimental results [ Fig. 4(c)]. The qualitative features of our calculations are not strongly sensitive to the value of J nn , provided that the interactions between pseudo-spins remain frustrated. Nevertheless, J nn = −0.64(4) K yields optimal fits, and we therefore use this value for calculations in the rest of this paper. We now consider the low-temperature state below T * . In this regime, most of the magnetic neutron-scattering intensity remains diffuse; however, weak magnetic Bragg peaks also appear on top of the magnetic diffuse scattering [ Fig. 4(b)]. The wavevector dependence of the scattering closely resembles observations for the Dy 3 Mg 2 Sb 3 O 14 , suggesting that a similar spin fragmentation process occurs in Ho 3 Mg 2 Sb 3 O 14 .
In particular, the divergence-full part of a spin-fragmented state describes an "all in, all out" (AIAO) long-range magnetic ordering involving only a small fraction ( µ/3) of the total magnetic moment [ Fig. 1(a)]. Rietveld refinements to the weak magnetic Bragg component of our T = 0.12 K data are in good agreement with the AIAO average magnetic structure in Ho 3 Mg 2 Sb 3 O 14 , consistent with a spin-fragmented state (see Appendix A). However, there are fundamental differences between our measurements and the CSF state observed in Dy 3 Mg 2 Sb 3 O 14 . First, the observed magnetic Bragg intensities are strongly reduced compared to the expected classical value. The magnetic Bragg intensity is proportional to the square of the ordered magnetic moment, which is µ/3 = 3.3 µ B per site for a CSF state in the absence of chemical disorder [ Fig. 1] [40,41]. In contrast, Rietveld refinements to our T = 0.12 K data indicate an ordered magnetic moment of only 1.70(3) µ B per Ho 3+ (see Appendix A). Importantly, Dy 3 Mg 2 Sb 3 O 14 has both a larger ordered moment (2.66(6) µ B per Dy 3+ at 0.20 K [41]) and a greater degree of site mixing than Ho 3 Mg 2 Sb 3 O 14 (see Section III). This suggests that chemical disorder cannot fully explain the observed reduction in ordered moment in Ho 3 Mg 2 Sb 3 O 14 , although it may be a contributing factor. Given the low temperature of our measurement (≈ 0.12 K), quantum fluctuations are the most likely alternative explanation. Second, the spin fragmented state in Ho 3 Mg 2 Sb 3 O 14 is accompanied by the persistence of continuous magnetic excitations down to at least 0.12 K. Furthermore, a distinct mode develops at ω ≈ 3 K when approaching T * and appears clearly separated from the elastic line in I(Q, ω) below T * [red arrows in Fig. 4(c)]. Above this mode, we observe a high-energy tail extending to ω ≈ 15 K that resembles the slow decay from the central peak in the paramagnetic phase. The presence of clear low-temperature spin dynamics over a wide energy range strongly contrasts with canonical classical Ising magnets such as Ho 2 Ti 2 O 7 [86] and Dy 3 Mg 2 Sb 3 O 14 [41], in which the spin dynamics is too slow to observe in neutron-scattering measurements at comparable temperatures [43]. The enhancement of inelastic scattering and reduction of the magnetic Bragg intensity thus provides experimental evidence for quantum excitations above a spin-fragmented ground state in Ho 3 Mg 2 Sb 3 O 14 .

VI. THEORETICAL MODELING
Our experimental results have revealed two key insights: that the spin Hamiltonian of Ho 3 Mg 2 Sb 3 O 14 realizes a frustrated transverse Ising model; and that a spin-fragmented state with quantum spin dynamics exists in Ho 3 Mg 2 Sb 3 O 14 at the lowest measurable temperatures. These results identify Ho 3 Mg 2 Sb 3 O 14 as an exotic frustrated magnet with unconventional properties. However, two important questions remain. First, how do frustrated interactions, transverse field, and hyperfine coupling conspire to generate the observed quantum spin dynamics? And, second, are these quantum dynamics primarily single-site fluctuations, or do they involve many-body quantum correlations as proposed for the QSF state in the introduction?
As a first step towards answering these questions, we perform a theoretical analysis of the spin Hamiltonian, Eq. (10). Because this Hamiltonian describes a many-body interacting quantum system with several interactions on the same energy scale, no single theoretical approach can provide a complete description of its properties. We therefore employ two complementary approaches. First, in Section VI A, we use exact diagonalization of small clusters to investigate quantum correlations in the putative QSF ground state. Second, in Section VI B, we develop a modification of Monte Carlo (MC) simulation in conjunction with real-space mean-field (MF) theory, which enables us to investigate finite-temperature properties on large clusters, with the compromise that quantum effects are now considered only at the single-site level. In both approaches, we identify how the hyperfine coupling competes with the transverse field to suppress quantum correlations.

A. Exact diagonalization
We use exact diagonalization (ED) of small clusters to explore the ground state of the model Hamiltonian, Eq. (10). The main results presented here are obtained from diagonalizing a 2 × 2 rhombus-shaped supercell with N = 12 sites; however, our conclusions are robust to the choice of different supercell geometries and the use of 18-site clusters (see Appendix E).
We start by considering a simple quantum model that contains the essence of Ho 3 Mg 2 Sb 3 O 14 : a dipolar kagome ice under transverse fields. In the absence of the hyperfine and exchange couplings, the system is controlled by a single tuning parameter h x /D. At the classical point h x /D = 0, we obtain 12 CSF states with energies that happen to be exactly degenerate for the 2 × 2 cluster. In a CSF state, the degeneracy grows exponentially with system size, and the ground state in the thermodynamic limit is long-range ordered with a √ 3 × √ 3 enlargement of the unit cell [39]. Our quantum model shows remarkably different behavior: an infinitesimal transverse field opens up a gap and stabilizes a unique ground state corresponding to a quantum superposition of all 12 CSF states with equal weights [ Fig. 5(a)]. The mixing of states out of the CSF manifold is negligible when h x /D is small, which can be seen by projecting the ground state wavefunction into the CSF subspace [ Fig. 5(b)]. To understand this state, we consider the pair correlation function of emergent magnetic charges on neighboring "up" triangles, Q Q m [ Fig. 5(b)], which would be unity for perfect spin fragmentation and zero for perfect paramagnetism. At small h x /D, we find that Q Q m is unity, suggesting that the divergence-full channel contains an ordered magnetic moment of µ/3, as in the CSF case [ Fig. 1(a)]. Moreover, our calculations show that the energy gap between the ground state and the first excited state scales as (h x /D) 6 until h x /D 0.2 [ Fig. 5(a)]. We take this as evidence that the leading-order quantum tunneling process appears at sixth order in perturbation theory and corresponds to flipping six spins around a hexagon, moving between degenerate CSF states [ Fig. 1(a)]. The existence of such lowenergy excited states hints at the presence of low-energy collective quantum dynamics in this model, although whether the excitation spectrum is gaped or gapless in the thermodynamic limit remains an open question. These results provide a theoretical foundation for the QSF phase by revealing that h x enables coherent quantum tunneling between classical kagome spin-ice states, yielding a quantum state that can be separated into an ordered divergence-full part and a quantum superposition of divergence-free channels.
As the transverse field is increased, the system undergoes a quantum phase transition at h x /D ≈ 0.2 into a nonperturbative regime, indicated by a peak in the second derivative of the ground-state energy [ Fig. 5(a)]. Notably, such a phase transition is absent in a transverse Ising model with only nearest-neighbor exchange interactions (see Fig. 9). In this phase, despite strong spin correlations, the energy gap significantly deviates from the (h x /D) 6 scaling of the QSF phase. The ground state quickly falls out of the CSF subspace, and consequently the emergent-charge correlation function shows a sharp drop. The system eventually crosses over into a paramagnetic phase in the limit h x D, corresponding to the single-ion limit.
For Ho 3 Mg 2 Sb 3 O 14 , there are three additional ingredients to consider: (i) coupling between kagome layers due to longrange dipolar interactions, (iii) nearest neighbor exchange interactions J nn , and (iii) hyperfine interactions A hf . While it is not possible to include inter-layer couplings because of the limited size of our ED calculation, we can nevertheless study the effects of J nn and A hf . In the classical case, J nn does not change the relative energy difference between different kagome spin-ice states. Therefore, a spin-fragmented state is still energetically favorable for ferromagnetic and weak antiferromagnetic J nn [39], as relevant here. To identify the effect of hyperfine coupling, we proceed as follows. First, we put Eq. (10) into a block-diagonal form where each block has a set of fixed nuclear spin numbers {I z 1 , I z 2 , · · · , I z N } and the hyperfine term is treated as a site-dependent local longitudinal field; this is possible because I z remains a good quantum number. Then, we find the global ground state of the system by diagonalizing the Hamiltonian in each of the 8 N nuclear spin blocks independently.
The interplay of hyperfine coupling with quantum correlations is shown in Fig. 5(c). Here, we use the parameters J nn = −0.64 K, h x = 0.85 K, and D = 1.29 K relevant for Ho 3 Mg 2 Sb 3 O 14 , and plot the nearest-neighbor spin and charge correlation functions as a function of the hyperfine coupling energy scale A hf I. In the absence of hyperfine coupling, the ground state is in the non-perturbative regime, with strong correlations between spins but weak correlations between emergent magnetic charges. For A hf = 0.319 K, as appropriate for Ho 3 Mg 2 Sb 3 O 14 , the ground state has I z i = ±7/2 for every site. With increasing A hf , the nearestneighbor quantum correlation-defined as the total correlation minus the classical correlation, σ z i σ z j − σ z i σ z jdiminishes quickly from 100% at A hf I = 0 to only 1% of the total correlation for the calculated ground state of Ho 3 Mg 2 Sb 3 O 14 (A hf I = 1.11 K). Simultaneously, the emergent charge correlation is greatly enhanced, with Q Q m = 0.71 for Ho 3 Mg 2 Sb 3 O 14 , corresponding to an ordered moment of 2.80 µ B per Ho 3+ [Fig. 5(c)]. These two observations together imply that quantum correlations in the groundstate of Ho 3 Mg 2 Sb 3 O 14 are strongly suppressed by hyperfine interactions. However, as the system is thermally excited to hyperfine levels with smaller |I z |, quantum correlations reappear: the most quantum state lives in the nuclear spin blocks where Fig. 5(c)]. Our ED results thus suggest that the low-temperature physics of Ho 3 Mg 2 Sb 3 O 14 can be understood in terms of a semiclassical ground-state resembling predictions from mean-field theory, yet with intrinsically quantum spin excitations.

B. Mean-field Monte Carlo simulations
Having established exact results for small spin clusters, we now develop a numerical method capable of approximating local quantum fluctuations in large spin configurations. The essence of our approach is to combine Monte Carlo simulations and exact diagonalization of the mean-field Hamiltonian at each site to obtain equilibrium configurations of static spins.
We begin by summarizing our approach, which we term "mean-field Monte Carlo" (MF-MC). The mean-field Hamiltonian of a given site i in a spin configuration is given by where is the longitudinal field that includes contributions from exchange, dipolar, and hyperfine interactions. Eq. (12) describes a 2 × 2 matrix in each sector n ∈ {1, · · · , 8} of the nuclear spin I z i . Diagonalizing these matrices yields the single-ion eigenstates in the same form as those shown in Fig. 2(d), i.e., |Ψ i n = ψ ± i , I z i . We obtain the static spin on site i as the diagonal matrix element of σ z i , For zero transverse field, the static spin σ z i is a classical Ising variable with unit magnitude on all sites. For nonzero transverse field, σ z i has a magnitude of less than unity, reflecting the nonzero probability of transverse-field-induced spin flips. The probability amplitudes of such spin flips are given by the off-diagonal matrix elements, To take account of this effect, in our Monte Carlo simulation we do not enforce σ z i = ±1; instead, we compute σ z i on-the-fly from mean-field states ψ ± i , I z i . Our protocol is as follows. We initialize σ z i with random uniformlydistributed values between −1 and 1. One MC step consists of selecting a random site i and proposing an update to the nuclear spin quantum number (I z i →Ĩ z i ) and electronic static spin, with the latter selected at random from three possibilities: (i) maintaining the static spin ( σ z i → σ z i ), (ii) flipping the static spin ( σ z i → − σ z i ), and (iii) going to one of the two new mean-field states obtained by diagonalizing the mean-field Hamiltonian, Eq. (12), for which both length and direction of the static spin may be changed ( σ z i → σ z i ). The move is accepted or rejected according to the Metropolis protocol (see Appendix F). Our method is identical to a classical MC simulation with single spin-flip moves in the limit of vanishing transverse field.
As in Section VI A, we begin by considering the simple case with A hf = 0. Fig. 6(a) shows a MF-MC phase diagram as a function of h x /D and J nn /D at low temperature (T = 0.12 K). We observe two trivial phases: first, for dominating h x , a paramagnetic ground state with σ z i = 0 is obtained as anticipated; and, second, for small h x and non-frustrated interactions (J nn /D 0), a conventional AIAO order is obtained, in which the mean field has the same magnitude on all sites. In contrast, frustrated interactions (J nn /D −1.5) favor non-trivial states in which the mean field-and hence the static spin-is spatially modulated [87]. Depending on the relative strength of h x compared to the pairwise interactions, this stabilizes two distinct phases. For h x D, we find a spin-liquid-like "one in, one out" (OIOO) phase characterized by a local constraint on every triangle: two out of the three spins are static with roughly equal moments and trace out closed loops in the kagome planes, while the third spin remains entirely dynamic due to the local cancellation of mean fields; hence, emergent magnetic charges are absent [ Fig. 6(b)]. Spin correlations of the OIOO phase gives rise to distinct star-like signature on Brillouin zone edges in momentum space, while magnetic Bragg peaks and pinch-point features associated with a spin-fragmented state are absent [ Fig. 6(c)]. For h x D, we find a phase resembling a CSF state but dressed with static moment modulation from site to site; we call this phase moment-modulated spin fragmented (MMSF). In this state, quantum fluctuations are manifest for each triangle in the form of one "long" static spin with a large magnitude, and two "short" static spins with smaller magnitudes [ Fig. 6(b)]. Remarkably, despite a broad distribution of the static spin length [ Fig. 6(e)], sharp magnetic Bragg peaks and diffuse scattering with pinch-point-like features coexist in momentum space [ Fig. 6(c)]. This observation provides a "smoking gun" for spin fragmentation, which can be understood in the following way. Because each triangle has a "one in, two out" or "two in, one out" arrangement of static spins, it supports well-defined emergent magnetic charges, Q k = ± i∈k σ z i , where i runs over three spins in a triangle k. Importantly, Q k form a staggered order with a small variation in magnitude [Figs. 6(b),(f)]. In the language of spin fragmentation, once we separate the divergence-full part of the MMSF phase (equivalent to the average value of Q k or the ordered moment of the AIAO structure), the remaining part is approximately divergence free and hence gives rise to diffuse scattering with pinch-point singularities.
We now compare our experimental data with MF-MC calculations using the parameters D = 1.29 K, h x = 0.85 K, and J nn = −0.64 K appropriate for Ho 3 Mg 2 Sb 3 O 14 . For A hf = 0, our MF-MC simulations predict a MMSF phase at low temperature [ Fig. 6(a)]. Setting A hf to the experimental value of 0.319 K enhances the classical spin correlations of the MMSF state, and increases the magnitudes of the static moments, ordered moments, and emergent magnetic charges [ Fig. 6(d-f)]. Our calculated heat capacity shows a significant improvement over the single-ion calculation and agrees qualitatively with the experimental data [ Fig. 3]. The calculation captures three major features: a broad shoulder between 1 and 10 K, a sharp phase transition into the spin-fragmented phase below 1 K, and a large nuclear heat capacity anomaly centered around 0.3 K. Our calculated diffuse-scattering pattern agrees well with our energy-integrated neutron scattering data at 4.2 K, yet discrepancies with experiment become increasingly significant as the simulation temperature is lowered [ Fig. 4(a)]. Most notably, at 0.12 K, the calculation is much sharper than the experimental data and strongly overestimates the intensities of the magnetic Bragg peaks. While the MF-MC calculated ordered moment of 2.85 µ B at 0.12 K is in good agreement with the ED result of 2.80 µ B , these values are considerably larger than the experimental result of 1.70(3) µ B [ Fig. 6(d)]. There are two likely explanations for this discrep-ancy. First, our models neglect the effect of Mg/Ho site disorder, but in fact approximately 3% of Ho sites are occupied by Mg in Ho 3 Mg 2 Sb 3 O 14 . Classical calculations suggest that this degree of site mixing can suppress the ordered moment by ∼ 0.5 µ B [41]; i.e., by a similar amount to quantum fluctuations. Second, the large specific-heat anomaly at low temperature makes thermalization of the sample challenging. According to our MF-MC calculations, if the true sample temperature is about 85% of T * , the ordered moment will be suppressed to the experimental value in the absence of site disorder.
Although our MF-MC calculations cannot describe collective excitations, the density of states of single-ion excitations shown in Fig. 4(c) provides insight into the dynamics of the MMSF state (see Appendix F). While such excitations are prohibited at low temperature in Dy 3 Mg 2 Sb 3 O 14 due to the vanishing off-diagonal neutron dipolar matrix elements, this is no longer the case in Ho 3 Mg 2 Sb 3 O 14 due to the presence of the transverse field [Eq. (15)]. In addition to a large central peak, our calculation shows two small inelastic peaks centered around ω = 5 K and 14 K that correspond to flipping short and long static spins, respectively. The relative amount of spectral weight in the elastic vs. inelastic channels is qualitatively consistent with the experimental data, but the distinct shoulder in the data appears roughly 3 K lower in energy transfer than in the calculation, and the calculation does not reproduce the continuous nature of the excitations. These observations may be a consequence of quantum correlations for the excited states, which are suggested by the ED results [Section VI A] but not included in the MF-MC model; hence, quantum calculations of the excitation spectrum would be interesting to explore in future studies.

VII. DISCUSSION & CONCLUSIONS
Our study demonstrates that Ho 3 Mg 2 Sb 3 O 14 realizes a frustrated quantum Ising magnet on kagome lattice with three competing energy scales: pairwise interactions between electronic spins, quantum tunneling via an intrinsic transverse field, and hyperfine coupling between electronic and nuclear spins. Our experimental study uncovers overdamped paramagnetic spin dynamics at high temperature and a spinfragmented state with reduced ordered moment at low temperature. The key difference with the CSF phase reported in isostructural Dy 3 Mg 2 Sb 3 O 14 is the observable low-energy spin dynamics in Ho 3 Mg 2 Sb 3 O 14 . These dynamics are a consequence of the intrinsic transverse field emerging from the lower crystalographic symmetry of the Ho 3+ site in Ho 3 Mg 2 Sb 3 O 14 compared to pyrochlore spin ices. Crucially, the transverse field is homogeneous, in contrast to random fields induced by chemical disorder. Therefore, a wideranging implication of our study is that symmetry lowering need not be a complicating factor in condensed-matter systems, but can instead enable the observation of simple models of quantum frustration.
Our results motivate further theoretical investigations of the interplay between spin fragmentation and coherent quantum tunneling. Our ED calculations suggest that the trans-verse field can generate ring-flip tunneling processes that connect degenerate CSF configurations and may stabilize a QSF phase in the thermodynamic limit. The presence of long-range dipolar interactions makes our model fundamentally different from quantum kagome Ising models studied previously [11,36,37,68,69] for which low-temperature states obtained at small transverse fields are continuously connected to the high transverse-field paramagnetic state in the absence of an external magnetic field [11,68,69]. Quantum Monte Carlo calculations are an especially promising approach for future research on quantum spin fragmentation, due to the absence of a sign problem for the model proposed here [88].
Finally, our results highlight the important role hyperfine interactions can play in frustrated quantum magnets based on rare-earth ions. It is usually assumed that nuclear spins have a negligible effect on the behavior of electronic spins, except in the case of singlet electronic ground-states where the hyperfine interaction may induce a cooperative ordering of the combined electronic and nuclear spin system [76,89,90]. Our results provide a more subtle example compared to this simple picture due to the highly-frustrated nature of spin interactions in Ho 3 Mg 2 Sb 3 O 14 . Here, hyperfine interactions induce static moments at the single-ion level, and enhance classical correlations at the expense of quantum correlations; yet, this does not drive the system towards a conventional long-range magnetically ordered state at low temperature. Finally, since all the stable isotopes of trivalent non-Kramers rare-earth ions have non-zero nuclear spin quantum number, hyperfine interactions may also have significant effects in other frustrated non-Kramers magnets, such as the quantum spin-ice candidates Pr 2 Zr 2 O 7 [22,91] and Pr 2 Hf 2 O 7 [24]. sored by the U.S. DOE, Office of Basic Energy Sciences, Scientific User Facilities Division. Identification of commercial equipment does not imply recommendation or endorsement by NIST.

APPENDIX A: STRUCTURAL AND MAGNETIC MODELS
The structural model of the s.s. sample of Ho 3 Mg 2 Sb 3 O 14 was obtained by Rietveld co-refinements to 50 K neutrondiffraction data collected using the HB-2A diffractometer [ Fig. 7(a)] and 300 K laboratory X-ray diffraction data [ Fig. 7(b)]. Refined values of structural parameters, and selected bond lengths and angles, are given in Table I. The canting angle of the Ising axes with respect to the kagome plane is 22.28(2) • from the co-refinement.
The average magnetic structure at low temperature was obtained by Rietveld refinement to energy-integrated neutronscattering data collected on the DCS spectrometer on our . Rietveld refinements to diffraction data for the s.s. sample. Co-refinements of the crystal structure to neutron and X-ray diffraction data are shown in (a) and (b), respectively. No obvious difference was observed between the X-ray diffraction pattern of s.g. and s.s. samples (not shown here). Refinement of the average magnetic structure to low-temperature magnetic diffraction data is shown in (c). In all panels, experimental data are shown as red circles, Rietveld fits as black lines, and difference (data-fit) as blue lines. Inset: illustrations of the all-in/all-out (AIAO) average magnetic structure within a unit cell and a single kagome layer.  (2)  s.s. sample [ Fig. 7(c)]. We isolated the magnetic Bragg scattering below T * by taking the difference between data measured at 0.12 K, and 0.4 K. The average AIAO magnetic structure belongs to the same irreducible representation as in Dy 3 Mg 2 Sb 3 O 14 , described by Γ 3 in Kovalev's notation [92], which is consistent with a spin-fragmented state [41]. The refined ordered moment is 1.70(3) µ B per Ho 3+ with a spin canting angle of 24.9 • with respect to the kagome plane.

APPENDIX B: POINT-CHARGE CALCULATIONS
Due to the low point symmetry at the Ho 3+ site, as many as 15 Stevens operators are required to describe the crystal-field Hamiltonian of the system [93]. The number of observables from the inelastic neutron scattering measurements is less than the number of unknown parameters, making conventional fitting procedures impracticable. To circumvent this problem, we calculated the crystal-field levels and wavefunctions from an effective electrostatic model of point charges using the software package SIMPRE [94]. The model considers eight effective oxygen charges surrounding a Ho 3+ ion, consistent with Rietveld refinements to the powder diffraction data (see Appendix A). The model is then adjusted numerically to match the measured crystal-field spectrum [95]. A similar method using a point-charge model has recently been applied to study the crystal-field spectrum in isostructural tripod-kagome compounds with Nd 3+ and Pr 3+ ions [96]. For The matrix element α = 0|Ĵ z |1 = 7.79 gives rise to a total magnetic moment of magnitude 9.74 µ B according to Eq. (4).
With eigenvalues E i and eigenvectors |i diagonalized from the above Hamiltonian, the magnetization along the local Ising axis is  Fig. 8]. At low temperatures, such as between 5 K and 1.8 K, deviations from the calculation become appreciable due to the development of spin-spin correlations that are not included in this model.

APPENDIX D: PARAMAGNETIC EFFECTIVE-FIELD FITS
We used an effective-field approach to calculate the inelastic neutron-scattering pattern in the paramagnetic phase, based on the Onsager reaction-field approximation [84]. Full details of this method are given in Ref. 83. In this approximation, the inelastic scattering function is given by where n = 3 is the number of Ho 3+ ions in the primitive unit cell, ω is energy transfer in K, and the susceptibility for each normal mode µ is given in Eq. (11) [80,97]. The magnetic structure factor is given by where Q is the scattering vector, r i is the position of magnetic ion i in the primitive cell, and z ⊥ i is its local Ising axis projected perpendicular to Q. The eigenvectors U iµ and mode energies λ µ are given at each Q as the solutions of where the Fourier-transformed interaction J ij (Q) = R J ij (R) exp (iQ · R) includes nearest-neighbor exchange and long-range dipolar contributions, and R is the lattice vector connecting atoms i and j. The dipolar interaction was calculated using Ewald summation [85]. The Onsager reaction field λ is determined by enforcing the total-moment sum rule The scattering intensities were calculated as where Q 0 = 0.4Å −1 and Q 1 = 1.6Å −1 , and where ω = 20 K, angle brackets here denote numerical spherical averaging, f (Q) is the Ho 3+ magnetic form factor [98], µ = 10 µ B is the total magnetic moment per Ho 3+ , and the constant C = (γ n r e /2) 2 = 0.07265 barn. The integrals were performed numerically.

APPENDIX E: EXACT DIAGONALIZATION
Due to computational limitations, our ED calculations were restricted to one tripod kagome layer. Cluster sizes of N = 12 and 18 with different shapes were studied under periodic boundary conditions. Long-range dipolar interactions were summed over periodic copies of the cluster cells up to a distance of 500 r nn .
The existence of a quantum phase transition at small h x /D as well as the the (h x /D) 6 dependence of the exciton energy are independent of size or shape of the clusters [ Fig. 9]. When truncating the dipolar interaction at the first nearest neighbor, our dipolar TIM becomes an exchange TIM model that has been investigated previously [68,69]. Our results are consistent with these studies for which the ground states obtained at low field are continuously connected to the high-field paramagnetic state with a h x /D dependence of the exciton energy. These results thus demonstrate that (i) the long-range dipolar interaction and the transverse field are two essential elements to realize the QSF state, and (ii) the QSF state is intrinsic to our effective Hamiltonian provided that the transverse field remains a perturbation to the dipolar interactions.
When considering the hyperfine interactions, Eq. (10) can be put into a block-diagonal form because I z remains a good quantum number. Then each block has a set of fixed nuclear spin numbers {I z 1 , I z 2 , · · · , I z N } and the hyperfine term is treated as a site-dependent local longitudinal field. Therefore, we find the global ground state by diagonalizing the Hamiltonian in each of the 8 N nuclear spin blocks independently. As expected, the ground states of the Hamiltonian lives in the nuclear spin blocks where I z i = ±7/2 for every site. The nuclear spin configurations of these ground states match with those electronic spins of the CSF states and thus possesses with the same degeneracy. In other words, the degeneracy of a classical Ising system is resorted in a quantum system through the nuclear spin channel. Recall that the ground state in the absence of hyperfine interaction is a superposition of all CSF basis states. The key effect of the hyperfine interaction is to promote the probability weight of one of the basis states over the others, which means it is suppressing quantum effects and driving the system to the classical limit, as discussed in the main text.

APPENDIX F: MEAN-FIELD MONTE CARLO SIMULATIONS
The key idea of our mean-field Monte Carlo simulation is to generalize the local states from the eigenstates of σ z to those of the single-ion mean-field Hamiltonian. These two sets of eigenstates are identical at the limit of h x = 0. The primary effect of the transverse field in our current calculation is to introduce quantum fluctuation at the level of a single site.
We keep track of the expectation values of σ z i , σ x i and I z i at each site during the simulation. The energy of the system is computed by replacing the spin operators in the Hamiltonian The relative difference of total energy between the Ewald summation and direct summation as a function of the splitting parameter α and reciprocal space cutoff n k . The real space cutoff is taken to be rc = min(a, b, c)/2, where a, b and c are the dimension of 5 × 3 × 2 orthorhombic supercell. The real space cutoff in the direct sum is rc = 1000rnn. The calculation is averaged over 50 random spin configurations. We take the reciprocal space cutoff n k = 10 and the best splitting parameter α = 7.78/ min(a, b, c). The same procedure is performed to find the best splitting parameter α = 8.10/ min(a, b, c) for the 10 × 6 × 4 simulation box with the same reciprocal-space cutoff.
Eq. (10) by their corresponding expectation values. The updates of static spin σ z i described in the main text are accompanied by corresponding changes in σ x i ; e.g., σ x i becomes − σ x i if the static spin is flipped. When a new mean-field state is proposed, we update both σ x i and σ z i . We used an orthorhombic 5×3×2 supercell of the crystallographic unit cell (see Appendix A) containing N = 540 sites in six kagome layers to obtain the heat capacity and the phase diagram shown in Figs. 3 and 6, respectively. A 10 × 6 × 4 simulation box containing N = 4320 sites was used to compute the Fourier transform of the static spin correlation function shown in Fig. 6(c), and the powder neutron-scattering patterns shown in Fig. 4. The long-range dipolar interaction was treated by Ewald summation with tinfoil boundary condition at infinity [99,100]. The interaction matrix was computed only once at the beginning of the simulation using the formulae for non-cubic unit cells given in Ref. 101. Suitable Ewald parameters were chosen by comparing with the result from direct summmation [ Fig. 10]. The canted local z-axis of 22.28 • obtained from the structure refinement was implemented. A simulated annealing process was performed during simulations. We started with random spin configurations at 20 K and cooled the system to 0.1 K with an exponential rate of 0.95. At each temperature, we performed 10000 equilibration sweeps, then collected samples after every 30 sweeps. The heat capacity shown in Fig. 3 was averaged over 30000 samples. One sweep consisted of N MC steps, where N is the number of sites in the simulation box.
The calculated magnetic scattering shown in Fig. 4(b) was obtained as the sum of static diffuse I diff (Q), Bragg I Bragg (Q), and inelastic I inelastic (Q) contributions, minus the high-temperature paramagnetic I para (Q) contribution, where the Bragg and diffuse contributions are calculated following Ref. 41. The inelastic contribution was given by The Q-averaged inelastic spectrum I(ω) was calculated as where E i = 2 (h x i ) 2 + (h z i ) 2 is the energy of the single-spin excitation and W i = (h z i ) 2 /((h x i ) 2 +(h z i ) 2 ) is the elastic spectral weight. The green line shown in Fig. 4(c) was obtained by convoluting I(ω) with the experimental energy resolution at the elastic line (FWHM ≈ 1 K).