Persistent spin dynamics in the pressurized spin-liquid candidate YbMgGaO$_4$

Single-crystal x-ray diffraction, density-functional band-structure calculations, and muon spin relaxation ($\mu$SR) are used to probe pressure evolution of the triangular spin-liquid candidate YbMgGaO$_4$. The rhombohedral crystal structure is retained up to at least 10 GPa and shows a nearly uniform compression along both in-plane and out-of-plane directions, whereas local distortions caused by the random distribution of Mg$^{2+}$ and Ga$^{3+}$ remain mostly unchanged. The $\mu$SR data confirm persistent spin dynamics up to 2.6 GPa and down to 250 mK with no change in the muon relaxation rate. They also reveal the power-law behavior of the spin-spin autocorrelation function that, along with the remarkable insensitivity of the system to the applied pressure, implies the relevance of structural randomness to the dynamic (spin-liquid) state in YbMgGaO$_4$.

Introduction. Spin-liquid states in frustrated magnets are nowadays actively studied as hosts for unconventional excitations representing magnetic monopoles [1,2] and other exotic quasiparticles [3,4]. Whereas several realworld materials show promising spectral features, their interpretation highlights the importance of multiple interactions that occur in real chemical compounds [5] and of the subtle structural features that may affect magnetic behavior therein [6].
Here, we focus on the spin-liquid candidate YbMgGaO 4 [7,8] that recently evolved as a unique triangular antiferromagnet with the robust three-fold symmetry, persistent spin dynamics, and a broad continuum of (potentially fractionalized) magnetic excitations. The crystal structure features triangular layers of the pseudospin- 1 2 Yb 3+ ions that are well separated by slabs of non-magnetic Mg 2+ and Ga 3+ . Thermodynamic measurements [7,8] and muon spin relaxation (µSR) [9] suggest the absence of magnetic order down to at least 50 mK. Weak spin freezing may occur around 100 mK [10], although it involves only a tiny amount of the magnetic entropy [7] and leaves no signatures in µSR [9].
Magnetic excitations of YbMgGaO 4 reveal a broad continuum that can be interpreted in terms of gapless spinons [11][12][13][14] or as arising from nearest-neighbor resonating valence bonds [15], the latter suggestion being particularly interesting, as it makes a direct link to Anderson's pioneering work at the outset of the spin-liquid research [16]. Material aspects come into play, though, when one takes into account the random distribution of Mg 2+ and Ga 3+ between the Yb 3+ layers (Fig. 1a). This structural randomness has drastic effect on the crystal electric field (CEF) levels of Yb 3+ [17,18] and likely eliminates magnetic contribution to the thermal conductivity [19].
Several of the recent theories of YbMgGaO 4 revolve around the models with quenched disorder, suggesting valence-bond glass [20] or randomly directed spin stripes [21,22] as the possible ground state. On the other hand, diffuse scattering [23] and magnetic excitations [11,15] can be described without invoking the structural randomness, only a weak distribution of the g-factors and exchange parameters was proposed [24], whereas theory anticipates the spin-liquid formation in a properly tuned triangular antiferromagnet even in the absence of randomness effects [25][26][27]. This makes it difficult to decide whether the structural randomness is integral to the physics of YbMgGaO 4 , or acts merely as a minor complication on top of the "genuine" spin-liquid physics caused by competing exchange interactions. The problem is in fact more general, because any real material contains certain amount of structural disorder that may be instrumental in stabilizing the spin-liquid state. For example, the kagomé spin-liquid compound herbertsmithite, Cu 3 Zn(OH) 6 Cl 2 [28,29], is prone to the Cu/Zn site mixing [30], but it can be driven to a magnetically ordered state under hydrostatic pressure [31]. This observation confirms that the spin liquid is triggered by a competition of exchange couplings, which is strongly influenced by pressure, whereas the Cu/Zn site mixing should not be affected at all. The stoichiometric pyrochlore compound Yb 2 Ti 2 O 7 [32] also shows a pressureinduced magnetic order above 0.1 GPa [33]. On the other hand, the off-stoichiometric samples of Yb 2+x Ti 2 O 7 are by far less sensitive to the applied pressure and do not develop any magnetic order up to at least 2.4 GPa [33].
Pressure evolution of a spin-liquid state can be a useful test, which we employ here to diagnose the driving force for the spin-liquid formation in YbMgGaO 4 . We demonstrate persistent and nearly unchanged spin dynamics up to pressures as high as 2.6 GPa, suggesting that structural randomness plays a significant role in this material. We further argue that the dynamics itself is typical of a system influenced by this randomness.
Single-crystal XRD. Room-temperature XRD data were collected at the ID15B beamline of the ESRF between ambient pressure and 10 GPa. A diamond anvil cell loaded with He gas and a small single crystal from the batch reported in Ref. 8 were used for the experiment. The position of Yb was fixed at (0, 0, 0), with its thermal ellipsoid refined anisotropically to account for local displacements caused by the random distribution of Mg and Ga. Altogether, 9 structural parameters were refined from about 120 symmetry-independent reflections collected at each pressure point.
The YbMgGaO 4 structure comprises edge-sharing YbO 6 octahedra layers separated by slabs of trigonal bipyramids jointly occupied by Mg and Ga (Fig. 1a). The R3m symmetry keeps all Yb-O distances equal but allows a trigonal distortion of the YbO 6 octahedra with the O-Yb-O angle β deviating from 90 • . Incidentally, this angle is equal by symmetry to the Yb-O-Yb bridging angle α between the octahedra (Fig. 1b). Two structural parameters, one distance and one angle, are thus sufficient to characterize both the local geometry of Yb 3+ and the geometry of the nearest-neighbor exchange pathway.
No changes in the crystal symmetry were observed in our experiment. Pressure evolution of the a and c lattice parameters demonstrates similar compressibility along the different crystallographic directions (Fig. 2a). This can be illustrated by fitting pressure dependen- cies with the second-order Birch-Murnaghan equation of state [34] resulting in the bulk moduli B 0,a = 151(3) GPa and B 0,c = 130(2) GPa and the volume bulk modulus B 0,V = 142(2) GPa. YbMgGaO 4 is thus more compressible than Yb 2 Ti 2 O 7 (B 0 = 219 GPa [35]) and other rareearth pyrochlores that typically feature the bulk moduli in excess of 200 GPa [36].
The Yb-O distances shrink by about 0.6% at 2.6 GPa (the highest pressure of our µSR experiment, see below) and by 1.7% at 10 GPa (Fig. 2b). The α = β angle shows a weak downward trend only. In order to confirm it, we relaxed the experimental structures using densityfunctional (DFT) band-structure calculations [34]. As DFT can not treat the Mg/Ga disorder explicitly, ordered structural models have to be used. This approximation leads to a constant offset between the DFT results and experiment. Nevertheless, not only the qualitative trends but also the slope are well reproduced (Fig. 2). We thus confirm that the Yb-O-Yb angle decreases under pressure. Compared to the ambient-pressure value, it changes by 0.07 • at 2.6 GPa and by 0.2 • at 10 GPa.
Local structure. Moderate changes of the average crystal structure are accompanied by a strong elongation of the Yb thermal ellipsoid. Whereas the in-plane displacements characterized by U 11 are not affected by pressure, the out-of-plane displacement component increases by 70% (Fig. 1d). This out-of-plane displacement has been linked to the local distortions of the YbO 6 octahedra caused by the random (and, generally, asymmetric) distribution of the differently charged Mg 2+ and Ga 3+ ions around Yb 3+ [17]. At first glance, the increase in U 33 implies a strong enhancement of the structural randomness under pressure, but the actual situation is more complex.
The displacement parameter U 33 is in fact gauging the displacements of Yb atoms relative to each other (∆z Yb , the distribution of Yb electron density along the c axis) and not relative to the oxygen cage. To confirm this, we performed DFT-based relaxations [34] for three possible scenarios, where the Yb 3+ ions are sandwiched between: i) two Ga 3+ ions on one side and two Mg 2+ ions on the other side (octahedra A); ii) one Ga 3+ and one Mg 2+ ion on each side (octahedra B); iii) two Ga 3+ ions on one side and a combination of Ga 3+ and Mg 2+ on the other side (octahedra C). Pressure has a strong effect on the Yb position in the octahedra A and nearly no effect on the Yb atoms in B and C (Fig. 1e). Surprisingly, this effect is a shift of the whole octahedron A along the c direction without any change in the octahedron itself. Using the geometrical parameters shown in Fig. 1c, we determine that for the octahedra A the deformation expressed by (d 2 − d 1 ) changes from 1.65 % at 0 GPa to 1.66 % at 10 GPa, whereas (β 2 − β 1 ) evolves from 4.29 • to 4.32 • [34]. The changes in the octahedra B and C are equally small.
The change in the position of the octahedra A is likely related to the accumulation of different charges above (Ga 3+ ) and below (Mg 2+ ) the Yb layer. But the key result at this juncture is that such a change reflects some re-arrangements within the Mg/Ga layers and does not affect any of the local physics. The YbO 6 octahedra undergo uniform compression and simply keep the deformation that they had at ambient pressure. The effects of CEF randomness should then remain unchanged.
µSR results. We now probe the evolution of spin dynamics. To this end, we use µSR, which is sensitive to internal fields as low as 0.1 G and probes the magnetic behavior locally. The experiments were performed at the GPD and Dolly spectrometers at the PSI at ambient pressure and at the GPD spectrometer under hydrostatic pressure [37] down to 250 mK on a polycrystalline sample. The data were collected in zero field (ZF) and in the longitudinal-field (LF) mode, where the magnetic field was applied parallel to the spin of the implanted muons. A double-walled MP35 pressure cell was used, similar to the design reported in Ref. 38, with Daphne oil 7373 as pressure-transmitting medium. The pressure value was determined by measuring the superconducting transition of a small piece of indium positioned next to the sample inside the pressure cell. Figure 3 shows the ZF µSR time spectra at ambient pressure and also at 2.6 GPa measured at 0.25 K and 15 K. The similar behavior of the µSR time spectra indicates no change of the ground state, and the absence of oscillations excludes pressure-induced magnetic order- ing in YbMgGaO 4 up to at least 2.6 GPa. To estimate the temperature dependence of the relaxation rate λ, we fitted the ZF-µSR time spectra by where f PC is the fraction of signal coming from the pressure cell, described by the A PC function from Ref. 37. The ensuing muon relaxation rate λ is temperatureindependent between 4 and 40 K (Fig. 3). Below 4 K, the increase in λ indicates the onset of spin-spin correlations that fully develop around 0.8 K, where λ flattens out and remains temperature-independent upon further cooling. This temperature evolution is essentially similar to the ambient-pressure µSR result reported in Ref. 9 and hardly changes at 2.6 GPa (Fig. 3). Even absolute values of λ are the same as at ambient pressure within the error bar.
Above 40 K, λ shows a steady decreasing trend described by an activated behavior λ = A + B e −∆/T with ∆ = 320 ± 20 K, which is reminiscent of the lowest CEF excitation energy of about 450 K [17]. This observation suggests that at high temperatures the relaxation is governed by an Orbach process [39] involving the excited CEF doublets of Yb 3+ . We also performed µSR in the presence of longitudinal field. Even a longitudinal field, which is more than 10 times higher than the local static field estimated from the low-temperature value of λ, is unable to decouple the muon relaxation, suggesting that the spins are dynamic in nature as has been seen at ambient pressure [9].
The LF spectra also show a universal scaling over three orders of magnitude [34] when plotted against t/H γ with γ = 0.35 at ambient pressure and γ = 0.85 at 2.6 GPa (Fig. 4). This scaling indicates a power-law behavior q(t) ∼ t γ−1 of the dynamic spin-spin autocorrelation function q(t) = S i (t) · S i (0) , as in spin glasses above the freezing temperature [40,41]. While this observation would not be inconsistent with the weak spin freezing detected below 100 mK in ac-susceptibility measurements [10], we note that neither ambient-pressure µSR [9] nor dc-magnetization [42] detect any spin freezing at this temperature. In fact, the power-law behavior of q(t) is not inextricably linked to the spin freezing. It can arise from critical fluctuations in systems influenced by structural disorder, such as intermetallics with the non-Fermiliquid behavior [43], but also in disorder-free systems in the vicinity of a quantum critical point [44,45].
Discussion. The absence of pressure effect on the dynamic spin state in YbMgGaO 4 could imply that pressure-induced structural changes are too weak to affect the magnetic couplings in this material, or that the spin-liquid formation is reinforced by the structural randomness, which itself is not affected by pressure. In the following, we argue that the structural changes in YbMgGaO 4 are in fact non-negligible, and the second scenario applies.
The changes in the Yb-O-Yb angle between 0 and 2.6 GPa may look minor (Fig. 2c), but they are not unimportant, because nearest-neighbor magnetic interaction is very sensitive to this angle, especially in the 90 − 100 • range where YbMgGaO 4 lies. For example, the spin-dimer compound TlCuCl 3 with the Cu-Cl-Cu angle of 95.6 • is in the gapped singlet state at ambient pressure but reveals pressure-induced magnetic order around 0.1 GPa [46][47][48]. In the Shastry-Sutherland magnet SrCu 2 (BO 3 ) 2 with the Cu-O-Cu angle of about 98.0 • , the pressure of 2 GPa reduces the spin gap, and drastically modifies magnetic excitations [49].
Anisotropic couplings in the Yb 3+ compounds should be influenced by the bridging angle too [50].
In NaYbO 2 , the smaller Yb-O-Yb angle leads to the twofold increase in the nearest-neighbor coupling [51,52]. Yb 2 Ti 2 O 7 shows pressure-induced magnetic order above 0.1 GPa [33], suggesting that even such a low pressure may significantly affect the couplings between the Yb 3+ ions. We estimated [34] that in Yb 2 Ti 2 O 7 the average Yb-O-Yb angle changes from 101.61 • at ambient pressure to about 101.67 • at 1 GPa. This change is comparable in magnitude to the 0.07 • change observed in YbMgGaO 4 between 0 and 2.6 GPa (Fig. 2c), suggesting that pressure effect on the couplings is not insignificant.
In the absence of structural randomness, the observation of the t/H γ scaling in the LF-µSR data would imply strong critical fluctuations in the vicinity of a quantum critical point and, consequently, high sensitivity of the system to the applied pressure, which is not observed. The only remaining option is that the scaling reflects (and confirms) the influence of structural disorder on spin dynamics. On the positive side, spin dynamics implied by the power-law scaling of q(t) is collective [40]. The continuous nature of magnetic excitations in YbMgGaO 4 [11,15,18] is then not a mere effect of structural disorder, but a signature of collective behavior and spin entanglement, as in "genuine" spin liquids, but, possibly, only on the length scale dictated by the structural randomness.
Conclusions. Hydrostatic pressure leads to a uniform compression of the YbMgGaO 4 structure with the reduction in the Yb-O-Yb angles, whereas local distortions of the YbO 6 octahedra are nearly unchanged. Spin dynamics is not affected by pressure and appears to be collective yet influenced by the structural randomness. This puts YbMgGaO 4 into the group of materials, where randomness effects are integral to the spin-liquid formation, but collective spin dynamics can be nevertheless observed.
We acknowledge ESRF for provision of beamtime. IEC thanks J. Jacobs for the He gas load. The work in Augsburg was supported by the German Science Foundation under TRR80 and by the Federal Ministry for Education and Research through the Sofja Kovalevkaya Award of Alexander von Humboldt Foundation (AAT).

Supplemental Material
Persistent spin dynamics in the pressurized spin-liquid candidate YbMgGaO 4

SINGLE-CRYSTAL XRD
High-pressure single-crystal X-ray diffraction on YbMgGaO 4 was measured at the ID15B beamline of the European Synchrotron Radiation Facility, Grenoble up to 10 GPa using monochromatic X-ray radiation (λ = 0.41114Å). Membrane-driven LeToullec-type diamond anvil cells (DACs) were used, equipped with Boehler-Almax anvils. Stainless steel was used as the gasket material, and helium was loaded as the pressuretransmitting medium. Diffraction patterns were collected with a Mar555 flat-panel detector using steps of 0.5 • oscillations over a total ω scan range of 76 • about the vertical axis. The pressures were measured using the ruby fluorescence method. Lattice parameter determination and integration of the reflection intensities were performed using the CrysAlisPro software [1]. Structures were refined using ShelxL [2] within the ShelXle [3] graphical interface.
Details of the data collection and structure refinement at 10 GPa are listed in Table S1. The refinements at other pressures were similar and can be found in the attached cif-files. We note that the ambient-pressure experiment was performed before loading the cell with the He gas. This may be the reason for the abrupt change in the Yb displacement parameters U 11 and U 33 between 0 and 0.78 GPa (Fig. 1d).
Pressure dependence of the lattice parameters and unit cell volume was fitted with the second-order Birch-Murnaghan equation of state, yielding B 0 = 142(2) GPa and V 0 = 253.12(7)Å 3 . Same equation was applied to the lattice parameters with the a and c values cubed. The EosFit routine [4] was used for the fitting.

General methodology
Structural information from single-crystal XRD can be extended by DFT-based structure relaxations. To this end, we performed non-spin-polarized calculations in the FPLO [5] and VASP [6] codes, and optimized atomic positions until the energy minimum was reached. The typical k-mesh included 64 points within the first Brillouin zone. Residual forces in the relaxed structures were below TABLE S1. Details of the single-crystal XRD data collection and structure refinement at 10 GPa.
3.34272(13) c (Å) 24.6033 (7)  The results of the relaxations are somewhat dependent on the methodology. In Table S2, we compare the influence of three factors: i) band-structure code; ii) exchange-correlation potential; and iii) treatment of the Yb 4f shell. The calculations are performed for the simplest ordered model of YbMgGaO 4 , where we restricted ourselves to the primitive cell of the R3m structure and reduced the symmetry to R3m, thus placing each Yb layer between the layers of Mg and Ga, respectively. This leads to a deformation of the YbO 6 octahedra that feature two Yb-O distances, d 1 = d 2 , and two O-Yb-O angles, β 1 = β 2 .
We conclude that the choice of the band-structure code affects the results of the relaxation. The calculations in VASP predict the too long Yb-O distances and, consequently, the too low angles. On the other hand, the relaxations within FPLO produce, irrespective of the methodology, the average distance (d 1 + d 2 )/2 and the average angle (β 1 + β 2 )/2 within, respectively, 0.01Å and 0.5 • from the experimental values. The larger deviations of the VASP results are probably related to the lower accuracy of the default pseudopotential for the Yb atoms, whereas FPLO does not rely on pseudopotentials and introduces no approximations to the crystal potential.
We note that all choices mentioned here (and even the VASP calculations) produce qualitatively similar results as function of pressure. The trends obtained from DFT are thus robust. In the calculations described below, we used FPLO, the GGA functional, and placed the Yb 4f states into the core, because it largely improves the convergence when bigger unit cells are considered.  (13) FIG. S1. Ordered structural models of YbMgGaO4 considered in this study. The structure shown in the left panel features weakly distorted YbO6 octahedra. One of them, the octahedron B (MgGaYbGaMg), is used to estimate relative displacements of the Yb atoms in Fig. 1e of the main text. The structure in the right panel covers the YbO6 octahedra with large deformations caused by the asymmetric environment, as in the octahedra A (MgMgYbGaGa) and C (GaGaYbMgGa). In both cases, the octahedra at z = 2 3 are only weakly distorted and can be used as reference for ∆z Yb , see text for details. VESTA software was used for crystal structure visualization [9].

Local structures
One difficulty in the ab initio modeling of YbMgGaO 4 is the random distribution of Mg and Ga, which cannot be included in the DFT calculation explicitly. We followed the approach of Ref. 10, where several ordered structures with the three-fold symmetry and YbMgGaO 4 composition are used to generate a gamut of representative local configurations of Yb 3+ . Those configurations were shown [10] to describe the broadening of crystal-field excitations and thus serve as a reasonable approximation to the real structure of YbMgGaO 4 .
In this work, we chose two of such ordered structures (Fig. S1) that correspond to the third and fourth structures from Fig. S5  the P 3m1 symmetry, with three nonequivalent Yb sites in each structure. In the left structure of Fig. S1, all Yb atoms develop a nearly symmetric local environment. We further used this structure to evaluate pressure dependence of the Yb-O distances and Yb-O-Yb angles shown in Fig. 2 of the main text. In this case, the octahedra are still deformed, but d 2 − d 1 < 0.005Å and β 2 − β 1 < 0.2 • indicate a negligibly small deformation. The "symmetric" structure (left panel of Fig. S1) is also used to obtain the displacements (∆z Yb ) of Yb within the octahedra B (GaMgYbMgGa).
The right structure of Fig. S1 covers the opposite scenario of the Yb atoms with a highly asymmetric local environment. It contains the octahedra A (GaGaYb-MgMg) and C (GaGaYbMgGa) showing the largest values of d 2 − d 1 and β 2 − β 1 . As explained in the main text, thermal displacement parameter of the average structure, U 33 , does not reflect deformations of individual octahedra and shows instead relative displacements of the Yb atoms with respect to each other. In the ideal structure, the Yb atoms should be at z = 0, 1 3 , 2 3 . To assess the effect of pressure on U 33 , we compare the positions of the Yb atoms within the octahedra A and C to the remaining Yb atom at z 2 3 that shows a rather symmetric local The fact that the octahedra in the reference layer are weakly distorted may look counter-intuitive, because their local environment (MgGaYbMgMg) is asymmetric too. This seeming discrepancy is probably caused by the simultaneous relaxation of the different octahedra within the same structure. Whereas the different octahedra could be explored independently by constructing larger supercells, this approach is unlikely to be fruitful, because in real structures different types of the octahedra coexist. Our model structures shown in Fig. S1 have an advantage that they describe not only the relaxation within the individual octahedra but also this coexistence effect.
Local deformations of the octahedra A, B, and C as function of pressure are shown in Fig. S2. We find that those octahedra undergo a nearly uniform compression and basically retain the deformation that they had at ambient pressure. This result implies that pressure has no visible effect on the structural randomness, even if the experimental U 33 parameter increases. Our model structures reproduce this increase in terms of ∆z Yb (Fig. 1e), but they do not show any significant changes in the local geometry apart from the gradual shortening of all Yb-O distances and an overall reduction in the β-angles.

Yb2Ti2O7
The calculations for Yb 2 Ti 2 O 7 followed similar methodology. According to Ref. 12, the low-temperature lattice parameter of the stoichiometric compound is very close to a = 10.0Å. Therefore, we performed calculations for this lattice parameter as well as for 9.9Å, and 9.8Å that, according to the Birch-Murnaghan equation of state with B 0 = 219 GPa and B 0 = 3.2 [11], correspond to pressures of 2.2 GPa and 4.4 GPa. The results shown in Table S3 suggest a consistent trend with dᾱ/dp = 0.055 • /GPa. Therefore, the change in the bridging angle for YbMgGaO 4 by 0.07 • between 0 and 2.6 GPa should be achieved in Yb 2 Ti 2 O 7 slightly above 1 GPa, whereas the pressure-induced magnetic order sets in already at 0.1 GPa [13] and requires a much smaller change in the angle.
µSR Figures S3 and S4 show the longitudinal-field (LF) decoupling experiments performed at ambient pressure and at 2.6 GPa. The universal scaling shown in Fig. 4 of the manuscript was verified using two different methods, similar to Ref. 14. First, the data points from all fields up to t = 6 µs [15] were arranged with increasing t/H γ for every value of γ. An empirical data mismatch function was calculated by taking the difference between the neighboring points and weighing them by the corresponding error bars. This mismatch function is defined as where N is the number of data points, and A i and δ i correspond to the asymmetry and error bar of the i-th data point, respectively. The lowest value of the mismatch function is obtained at γ = 0.35 (0 GPa) and γ = 0.85 (2.6 GPa) that produce the universal scaling over at least three orders of magnitude in t/H γ . We note that such a scaling may not hold in the whole time and field range -e.g., at very short times comparable to the width of the muon pulse [14] or in very high fields that affect spin dynamics -but the scaling over three orders of magnitude is comparable to all earlier observations [16,17] and serves as a robust evidence of the power-law behavior of the spin-spin autocorrelation function that, in turn, indicates collective spin dynamics. Additionally, we checked the scaling behavior by fitting individual relaxations with a stretched exponential, e −[λ(H)t] β , supplied with the time-independent background. This background originates from the thin copper plate sample holder at ambient pressure in the Dolly Heliox cryostat and from the pressure cell contribution at 2.6 GPa above 100 G. The values of β = 1 and β = 0.6 were used for the ambient-pressure data and 2.6 GPa data, respectively. The λ(H) then shows the power-law behavior, λ(H) ∼ H −γ , with γ 0.3 (0 GPa) and 0.8 (2.6 GPa) in good agreement with the direct scaling of the LF data. * altsirlin@gmail.com