Giant thermal magnetoconductivity in CrCl$_3$ and a general model for spin-phonon scattering

Insulating quantum magnets lie at the forefront both of fundamental research into quantum matter and of technological exploitation in the increasingly applied field of spintronics. In this context, the magnetic thermal transport is a particularly sensitive probe of the elementary spin and exotic topological excitations in unconventional magnetic insulators. However, magnetic contributions to heat conduction are invariably intertwined with lattice contributions, and thus the issue of spin-phonon coupling in determining the spin and thermal transport properties becomes more important with emergent topological magnetic system. Here we report the observation of an anomalously strong enhancement of the thermal conductivity, occurring at all relevant temperatures, in the layered honeycomb material CrCl$_3$ in the presence of an applied magnetic field. Away from the magnetically ordered phase at low temperatures and small fields, there is no coherent spin contribution to the heat conduction, and hence the effect must be caused by a strong suppression of the phonon thermal conductivity due to magnetic fluctuations, which are in turn suppressed by the field. We build an empirical model for the thermal conductivity of CrCl$_3$ within a formalism assuming an independently determined number of spin-flip processes and an efficiency of the phonon scattering events they mediate. By extracting the intrinsic phonon thermal conductivity we obtain a quantitative description at all fields and temperatures and demonstrate that the scattering efficiency is entirely independent of the field. In this way we use CrCl$_3$ as a model system to understand the interactions between spin and phonon excitations in the context of thermal transport. We anticipate that the completely general framework we introduce will have broad implications for the interpretation of transport phenomena in magnetic quantum materials.


I. INTRODUCTION
Transport measurements form one of the three pillars of experimental condensed matter physics. In insulating crystalline systems, the thermal conductivity κ (T ) ranks as one of the most valuable probes for investigating the low-energy excitations [1]. Unlike thermal equilibrium quantities such as the specific heat, c(T ), κ (T ) is a steady-state transport property and thus contains fundamental information about the itinerant characteristics of a system, most notably the relaxation times and scattering strengths of the low-energy excitations. In the field of insulating quantum magnets, lowdimensional spin systems may host a wide range of exotic ground states and κ (T ) has long been one of the most important probes of their unconventional spin excitations [2][3][4][5][6][7][8][9][10]. Because the lattice phonons invariably constitute a strong and relatively well-characterized contribution to thermal transport, κ (T ) measurements offer particular insight into the dominant spin-phonon scattering mechanisms. Even in materials whose thermal conductivity is phonon-dominated, meaning that there is no significant heat transport due to coherent magnetic modes, a strong field dependence of κ (T ) may still be present due to destructive effects of the spin sector on the phonon transport.
CrCl 3 is an insulating, layered, honeycomb-lattice compound and has attracted considerable recent attention from two independent lines of research. One concerns the "candidate Kitaev" material α-RuCl 3 [11][12][13][14][15], whose proximity to Kitaev physics may be gauged from the nature of its magnetic excitations [16][17][18][19][20]. While κ (T ) measurements have not been able to provide conclusive evidence of fractionalized fermionic spin modes in α-RuCl 3 , one of the primary reasons why the issue remains open concerns the role and indeed the nature of spin-phonon scattering, which is manifestly strong over a wide range of temperatures [10,21]. CrCl 3 is the 3d transition-metal structural analog of α-RuCl 3 , and as such represents the latter material in the absence of significant spin-orbit coupling. Although one may fear that this removal of Kitaev character removes any connection between the two systems, we will show here that CrCl 3 presents a test case for spin-phonon scattering effects that are at least as strong as any comparable phenomena in α-RuCl 3 .
The second avenue leading to CrCl 3 as a key material to understand is its position in the structural series CrX 3 , where X = Cl, Br, I is a halide. Structurally, the chromium trihalides are van der Waals materials, allowing them to be cleaved easily and prepared in mono-or few-layer forms that show strong differences in their physical properties. Magnetically, the honeycomb layers have ferromagnetic (FM) in-plane interactions and indeed both CrI 3 and CrBr 3 are bulk ferromagnets. This magnetic character has therefore promoted a keen interest in spin and lattice control in the context of topological magnonics [22], spintronics, and magnetoelectronics [23,24]. CrCl 3 is a historical "mixed FM/AF" system, with antiferromagnetic (AF) interlayer interactions ensuring an AF ground state of antialigned FM layers [25][26][27][28][29], and because of the rather low field scale (μ 0 H s 2 T, independent of direction) for complete spin polarization, is an excellent candidate for studying magnetic-field effects on spin and lattice transport.
In Fig. 1, we show our key experimental result, that the in-plane thermal conductivity of CrCl 3 is greatly enhanced by a magnetic field. We comment for clarity that by "inplane" we refer to measurements of κ (T ) performed with the temperature gradient ∇T , oriented in the ab plane and with the field also applied in the plane (H ∇T ). Qualitatively, this field dependence of κ (T ) resembles the behavior typical of magnetoresistance in magnetic conductors, where the electrical conductivity often increases as the field is increased [30]; this effect is caused by the field-induced reduction of the spin-dependent scattering and is often most prominent in the vicinity of the magnetic ordering transition. In CrCl 3 , the origin of this giant thermal magnetoresistance lies in the very strong field-induced suppression of the spin-phonon scattering, and one may also observe in Fig. 1 that it is particularly prominent at the ordering transition (T N ). Quantitatively, the phenomenon is anomalously large and is particularly unusual in that it extends over essentially the entire range of temperatures shown in Fig. 1. Here we will show that it can be modelled accurately with a minimal number of assumptions and empirical parameters.
In many cases, the phonon thermal conductivity is well captured by the highly refined Debye-Callaway (DC) model [31,32], in which all of the most important mechanisms for phonon scattering are considered over a wide range of temperatures. One obtains the form where ω is the phonon frequency, D is the Debye temperature, v s is a characteristic average phonon velocity, and the integration variable is x =hω/k B T . In the relaxationtime approximation, it is assumed that all possible scattering mechanisms contribute independently to the phonon scattering time, τ (ω, T ), and hence where the four relaxation times account respectively for boundary scattering, point-defect scattering, umklapp scattering, and resonant scattering due to impurities, magnetic excitations, or other collective modes. A similar formalism can be adopted for the direct contributions to κ (T ) of welldefined magnon modes. This type of model has been used to obtain a quantitative account of the thermal conductivity in a number of low-dimensional spin systems where both the spin excitations and their phonon-scattering effects, appearing in the τ −1 res term, can be characterized accurately [3,5,6,33]. Despite the experimentally verified success of the DC model for certain cases, in practical applications the model of Eq. (2) requires at minimum seven to ten fitting parameters to reproduce even rather smoothly varying κ (H, T ) curves. The microscopic implications of each term are often very difficult to verify independently, unless the respective fitting parameters can be compared among closely related materials (for example, by doping or elemental substitution). Most importantly for the CrCl 3 problem, there are no well-defined spin excitations over most of the (H, T ) phase diagram, and hence no possibility of describing the giant spin-phonon scattering within a τ −1 res term. Clearly the form of this term in Eq. (2) is too restrictive to capture the rich spectrum of possible interactions between the phonons and magnetic excitations of a quantum magnet, particularly if the latter are fractionalized. Thus we will introduce a more general approach to modeling the thermal conductivity of a magnetic insulator, a task in which we will be phenomenological but quantitative.
The focus of our contribution is to describe the heat conduction of magnetic insulators in the common situation where this is governed mostly by phonons, but subject to a spin-phonon scattering to which multiple mechanisms may contribute. In such a case, it would be highly desirable to have a method of understanding the magnetic scattering of phonons without invoking either the microscopic details of the scattering mechanism or system-specific characteristics such as the phonon and magnon (or spinon) dispersion relations. Here we present a phenomenological model to quantify the T and H dependence of κ, for the worked example of CrCl 3 , by considering the phonon heat conduction in the presence of scattering by magnetic degrees of freedom. The parameters determined empirically in our model enable us to quantify the dominant scattering mechanisms regardless of the energy scales of the phonons and of the magnetic excitations.
The structure of this article is as follows. In Sec. II, we summarize briefly our samples and experimental methods. In Sec. III, we show the results of our κ (H, T ) observations at all measured fields and compare these with measurements of the magnetization and specific heat. In Sec. IV, we present the details of our empirical modeling procedures for the physically different regimes and extract the quantities required to describe spin-phonon scattering in CrCl 3 . Section V provides a discussion and conclusions.

II. MATERIAL AND METHODS
Thin, purple, platelike single crystals of CrCl 3 were grown by chemical vapor transport [34,35]. CrCl 3 is known [36] to have a rhombohedral low-temperature crystal structure composed of hexagonal lattices of Cr 3+ ions in the ab plane, whose ABCĉ-axis stacking is ensured by only rather weak (van der Waals) structural interactions. Magnetically, as noted in Sec. I, the S = 3/2 (high-spin) Cr 3+ ions in the honeycomb layers have FM interactions of strength J = −5.25 K [28], which is a consequence of the near-90 • Cr-Cl-Cr geometry (edge-sharing CrCl 6 octahedra), while the interlayer interactions (J ) are AF.
The magnetization m(T ) and heat capacity c(T ) were measured over a range of temperatures from 2 to 100 K and of applied magnetic fields up to 18 T using, respectively, Quantum Design MPMS and PPMS systems. The in-plane longitudinal thermal conductivity, κ xx ≡ κ, was measured on an as-grown sample of dimensions 2×4×0.5 mm using a single-heater, two-thermometer configuration in steady-state operation with the field applied in the ab plane and in the direction of the thermal gradient (∇T H ∈ ab); limited κ measurements were also performed with the field normal to the plane. The difference in absolute temperatures across the sample was set never to exceed 5% of the bath temperature throughout the entire T range of the measurements. All thermometry was performed using Cernox resistors, which were precalibrated individually and in situ under the maximum applied fields of each instrument.

III. EXPERIMENTAL MEASUREMENTS
A. Thermal conductivity κ(T ) Figure 1 shows the complete picture of κ as a function of temperature for all values of the applied in-plane field, H, that we measured. At zero field (ZF), κ (T ) exhibits a somewhat flat maximum at 20-25 K with a gentle decline to higher temperatures; at lower T , it has a marked plateau-type region at and below the magnetic ordering temperature, T N = 14.2 K. As H is increased, it is clear that fields on the order of μ 0 H = 1 T have only a minor effect on κ (T ). However, beyond 1 T, the applied field causes a dramatic increase of κ (T ) at all temperatures and the development of a strong and sharp peak at T p 20 K; at 18 T, the maximum exceeds its ZF value by a factor of 2.2. As we will show below, in fact κ (T ) is almost completely independent of the direction of the applied field, indicating a minimal magnetic anisotropy. At high fields, the peak at T p and the line shape on both sides of it, which retains only a minor remnant of the plateau at T N , are characteristic of phonon-dominated thermal conductivity. The value of T p varies little as H is increased. Its physical origin lies in a crossover between the dominant phonon scattering mechanisms. At T > T p , umklapp processes dominate and τ −1 U is the largest term in Eq. (2); for T < T p , defect-(τ pd ) and boundary-scattering processes (τ b ) take over. It is reasonable to assume that κ (T ) at μ 0 H = 18 T is closest to reproducing the purely phononic response, κ ph (T ), of CrCl 3 , and we return to this topic in detail in Sec. IV C. As the inset of Fig. 1 makes clear, at this field the low-T κ (T ) exhibits a T 3 dependence, suggesting that the thermal conduction is due to ballistic transport of acoustic (linearly dispersive) phonons.
At all lower fields, including zero, κ (T ) reflects a systematic suppression due to additional spin-dependent phononscattering processes over essentially the entire T range. Only at very low temperatures (T < 5 K) does the ZF κ exceed that at all other fields, because this is where the contributions of coherent magnon excitations in the magnetically ordered phase become important. From the inset of Fig. 1, the thermal conductivity in this regime has no simply characterized form, and may be a consequence of comparable magnon and phonon contributions, at least one of which does not show ballistic transport [37]. By contrast, at all temperatures above 5 K the dominant effect of the spins is not an additive contribution, from three-dimensionally coherent excitations, but a destructive effect that we assume is due to scattering of the phonons by incoherent spin fluctuations. It is the field-induced suppression of these fluctuations that brings about the striking enhancement we observe in κ.

B. Magnetization
To understand this interplay of lattice and magnetic excitations, we first consider the magnetic properties of CrCl 3 . The magnetization m(H ) is shown in the inset of Fig. 6 for several different temperatures. As noted in Sec. I, a moderate field μ 0 H ≈ 2 T, applied in the plane, is sufficient at T < 5 K to rotate the FM planes against the AF interplane interaction and drive the system to saturation. We comment that the field scale for this process suggests, in contradiction to the conclusion of Ref. [28], that the interplane J is of order 1 K. We deduce a saturation moment, m s = m(T → 0), of 2.88 μ B per Cr 3+ ion. The same value of m s is obtained by applying the same field in the out-of-plane direction, demonstrating that in CrCl 3 , unlike RuCl 3 , the magnetic anisotropy is extremely weak [35]. This m s value is fully consistent with the spin contribution expected for a single Cr 3+ ion of S = 3/2 when the g factor appropriate for a magnetically isotropic system, g = 2, is assumed.
For perspective on the relation between κ (T ) and m(T ), in Fig. 2, we show κ (T ) on the same temperature axis as the field-normalized magnetization, m(T )/H ≡ χ (T ), for fields of 0 and 1 T in κ and 0.01 and 1 T in χ . This low-field regime is where the plateau in κ (T ) around T N is most pronounced, and the field is not sufficient to suppress the spin-induced scattering of the phonons conducting the heat. For very small fields, it is clear that χ (T ) becomes very large in the region of the ordering transition, indicating strong magnetic fluctuations and hence the origin of the κ plateau. This behavior is suppressed considerably at 1 T, where the plateau begins to lose its form. Above 25 K, the residual m scales precisely with H.

C. Magnetic specific heat
To obtain further insight into both the phonons and their spin-mediated scattering, next we consider the specific heat as a function of T and H. Figure 3(a) shows our measurements of the total specific heat, c(T ), which are fully consistent for all H with previously reported results [35]. To estimate the phonon contribution, c ph (T ), we employ a three-dimensional (3D) Debye-model fit of heat-capacity data for the nonmagnetic analog ScCl 3 [38]. The resulting c ph (T ), shown by the green line in Fig. 3(a), agrees well with our data for CrCl 3 in the high-temperature regime, which as the inset illustrates is T 100 K.
We attribute the remaining, strongly field-dependent portion of the specific heat to the magnetic degrees of freedom, meaning we define c mag (T ) = c(T ) − c ph (T ). Figure 3(b) shows c mag (T ) on a logarithmic temperature axis. At ZF, there are clearly two peaks, a sharp, λ-type anomaly at the Néel transition, T N = 14.2 K, and a rounder maximum at 17-18 K. When the applied magnetic field is increased to 1 T, the λ anomaly is completely lost, and in fact a detailed study of the low-field evolution of this feature [35] found that a field of 0.2 T is sufficient to suppress this indicator of the ordering transition. As the applied field is increased, the location T max (H ) of the broad maximum shows an abrupt initial shift towards higher temperatures, before increasing more slowly with μ 0 H to a value of 36 K at 14 T.
Quite generally, the broad maximum in c mag (T ) fingerprints the energy scale of the dominant local spin-flipping processes in the system. In CrCl 3 , this is the energy cost for reversing a single spin in the FM honeycomb layers (in the rhombohedral structure one would anticipate an energy of 3|J| + J ). The initial increase in T max (H ) can be understood as a straightforward reinforcement of this magnetic stiffness while the applied field competes with J to reorient the FM layers. The slower increase at μ 0 H > 2 T corresponds to a competition of the field with 3|J|, which in Ref. [35] was formulated as a progressive development of ferromagnetic correlations in the plane. Because of the low spin-coercivity in an applied field, FM alignment of very large domains becomes thermodynamically favorable at higher T . It is important to note that T N , the temperature for 3D magnetic order, almost coincides with the broad maximum at ZF. Thus, despite the quasi-2D nature of the structure of CrCl 3 , the AF state at ZF is magnetically 3D, with T N of the same order as the Curie-Weiss temperature (θ CW 30 K both from our data below the structural transition at 240 K and from that of Ref. [35]).
While one may fear that our estimated c mag (T ), as a small difference between two larger numbers, is subject to significant errors, a complete validation of its accuracy may be obtained by integration. The magnetic entropy, . H is applied in the ab plane in every panel, and also along the c axis at T = 13.5 and 21 K, which highlights the very low magnetic anisotropy. For all temperatures T > 4 K, the applied field causes a monotonic enhancement of κ. Only at T < 4 K is a different behavior observed as a consequence of the coherent magnonic contribution: κ is enhanced by a small applied field (which causes magnon stiffening) but suppressed by all higher fields (where the energy cost of magnon excitation becomes prohibitive).
Clearly at temperatures beyond 50 K, s mag (T ) for all fields displays a smooth and accurate approach to the limit k B ln 4 expected for the four-level system corresponding to free S = 3/2 spins. This result, which contrasts strongly with the entropies shown in Ref. [35], indicates that our approach of deducing the lattice contribution from the material ScCl 3 provides a quantitatively accurate estimate of c ph (T ).

D. Thermal conductivity κ(H )
The clearest way to gauge the effects of the magnetic state on the thermal conductivity is to consider the isothermal H dependence of κ. In Fig. 4, we show for a range of temperatures the fractional change κ/κ 0 = [κ (H ) − κ (0)]/κ (0) caused by an applied magnetic field. It is apparent immediately that the generic field response is a monotonic increase in κ/κ 0 , except at the lowest temperatures.
Although our T = 1.6 K data are something of an outlier for our present purpose, which is to discuss the generic behavior beyond the ordered regime, we focus first on the low-T case. Here κ (H ) increases initially, peaking around 2 T before decreasing to negative values of κ/κ 0 and becoming H-independent at fields μ 0 H > 8 T. This behavior, which is observed for T < 4 K, indicates a coherent magnonic contribution to heat conduction that is comparable to the phononic one. The coherent spin-sector contribution to thermal conductivity has been the focus of numerous studies [2,3,6-8,10] and can in general be described as an independent additive term beyond the phonon contribution, i.e. κ tot = κ ph + κ mag . At low fields, the initial rise of κ (H ) can be explained by the magnon stiffening that occurs when μ 0 H overcomes J to orient all of the layers ferromagnetically, which results in an increase of the magnon propagation speed and hence of κ mag [39]. However, this contribution decreases as soon as the field-induced magnon gap exceeds the measurement temperature, at which point the magnon population is suppressed exponentially, leaving a smaller and largely H-independent thermal conductivity, κ tot → κ ph . The coherent magnon contribution also diminishes quickly as T increases, and in fact κ/κ 0 changes character well below T N (e.g., T = 8.5 K in Fig. 4).
At all temperatures T > 4 K, it is clear that the leading role of magnetic fluctuations is as scattering centers, rather than as participants in any kind of coherent heat conduction. While this type of behavior has been observed in certain quantum magnetic materials at specific temperatures, where it can be described by a τ −1 res term in Eq. (2) [5,6,9], systems in which it appears across the full range of temperatures are not widely known. Nonetheless, κ increases monotonically with field in all the relevant panels of Fig. 4; we remark again that this behavior is largely independent of the field direction. While it remains the case that increasing H suppresses the spin fluctuations, its effect is a uniform suppression of the spin-phonon scattering by an effective reduction in the density of scattering centers, thus bringing the system closer to purely phononic heat conduction.

IV. EMPIRICAL MODEL FOR κ(H, T )
To model the thermal conductivity in the presence of such a large and obviously destructive spin-phonon interference effect, we base our analysis on the thermal conductivity, κ ph (T ), due to phonons. We first reexpress Eq. (2) in the form i.e. we assume that the effective phonon scattering rate can be separated into a part τ −1 0 ≡ τ −1 0 (T ) containing all the fieldindependent terms in Eq. (2) and a part τ −1 mag ≡ τ −1 mag (H, T ) containing all phonon scattering processes involving magnetic degrees of freedom.
To capture the behavior of κ (H, T ) phenomenologically, we assume that the relative scattering rate, τ −1 mag /τ −1 0 , will depend on the population of magnetic fluctuations, which we denote n mag , and on a dimensionless coupling constant describing the strength, or effectiveness, of the phonon scattering processes, λ(H, T ). Thus the key equation underpinning our empirical treatment is that the thermal resistivity will be given by where κ ph , the phonon thermal conductivity in the absence of spin-phonon interactions, is H-independent and would be recovered at high applied fields. By some straightforward algebra, one may verify that λn mag = τ −1 mag /τ −1 0 = τ 0 /τ mag . We normalize n mag to the total number of Cr 3+ spins so that it represents a fractional density, and thus the field dependence of κ is encoded in two unitless parameters, 0 < n mag < 1 and 0 < λ.
Focusing now on the spin fluctuations responsible for phonon scattering, because these change character across the (H, T ) phase diagram, we begin by subdividing this into the two regimes shown schematically in Fig. 5. Qualitatively, at high fields and low temperatures one expects the fluctuations to be well-defined but 2D spin-wave excitations within the FM layers, whereas for temperatures high relative to the field one expects random fluctuations of weakly coupled spins. Specifically, for CrCl 3 , we define region I as covering low temperatures (T < T N ) for fields H 2 T and T < 30 K at our highest measurement field; here we will find that n mag corresponds closely to the density of magnon excitations originating within the honeycomb layers, which are highly spin-polarized and strongly ferromagnetically correlated, and hence remain coherent 2D entities outside the AF phase. Region II covers temperatures T > T N for low fields, and at our higher fields temperatures T > 40 K; here n mag corresponds to the average density of free spins that are not aligned with the field direction, and as such is well described by a Weiss-field picture, within which the net magnetization is determined by conventional paramagnetic behavior. As Fig. 5 makes clear, the AF ordered state occupies a very small region of the (H, T ) phase diagram, and because the contributions of coherent 3D magnon excitations to κ are also small (Fig. 1), the ordered regime is not the focus of our study. However, we comment that short-ranged magnetic correlations in the vicinity of the ordered state (gray colors) will act to limit the applicability of the approximations we apply below.
By considering the magnetic scattering in the high-H/T I and low-H/T regions II of Fig. 5, we will minimize the ambiguity at intermediate H/T ratios. In both regions, the fractional density of magnetic fluctuations, either of 2D spinwave excitations or of freely fluctuating spins, can be fixed accurately to known limits. In region I, this limit is the fractional deviation, 1 − m(H, T )/m s , of the net magnetization from its high-H/T saturation value, which is m s ≡ m(H → ∞). In region II, it is the fractional polarization, m(H, T )/m s , which vanishes in the low-H/T limit.
It is apparent from Eq. (4) that the type of phenomenological approach we adopt requires a reliable estimate of the pure (field-independent) phonon thermal conductivity κ ph (T ) to be useful or even viable. If the dependence of κ on the applied field is sufficiently weak at strong fields, many authors [7,40] use the κ (T ) data at their highest available magnetic field as a measure of κ ph . In Sec. IV C, we will use our modeling procedure to obtain a field-independent κ ph (T ) by extrapolating from our κ (H, T ) measurements (shown in Fig. 1), and hence will determine how close our κ (μ 0 H = 18 T) data are to saturation.

A. Region I: dominant magnon scattering
Our magnetization data (inset Fig. 6) show that, for temperatures T < 30 K, fields μ 0 H 5 T are strong enough to polarize more than 50% of the Cr 3+ spins. This temperature and field range lie well within the "spin-flop" regime [28,35], where the spins are forced into a quasi-FM-ordered state of the honeycomb layers. To model this regime, we assume 2D FM spin-wave excitations with the field-gapped dispersion relationhω(k, H ) = 1 2J (ka) 2 + gμ B H, valid at small k, where a is the lattice constant andJ is an effective magnetic interaction strength. The population of magnetic excitations per unit volume in d dimensions may then be estimated using the expression [41] n mag (H, From the structure of this expression, only wave vectors k near the zone center contribute significantly to n mag , and thus we take the upper limit of the integral to infinity; it is also not necessary to specify the departure of ω(k, H ) from a simple quadratic at higher k values. α is a unitless constant set by the dimension, d, and the normalization. Li s (x) is the polylogarth- mic function of order s, which has an infinite series expression that is easy to evaluate numerically. As Fig. 6 makes clear, this approach provides an excellent account of the H and T dependence of m, effective even at temperatures (21 and 32 K) outside region I as it is represented in Fig. 5. In Fig. 6, the only fitted parameter is the effective interaction,J = 13.1 K; clearlỹ J ≈ 3J is again close to the characteristic energy scale of the 2D FM layers at ZF. Here we have treatedJ as a constant, anticipating that fields up to 7 T cause only weak changes in the effective in-plane stiffness, certainly when compared with their gapping effects (this is to be contrasted with the "magnon stiffening" discussion in Secs. III C and III D, where the field dominates J ).
At high H/T , when the total magnetization remains a significant fraction of m s , we equate the measured deviation directly with n mag , the population of 2D magnons obtained using Eq. (5) with d = 2. The agreement remains good for all temperatures T < 2T N when μ 0 H 2 T. Given the discrepancy between the in-plane FM interactions and the outof-plane AF ones, whose 3D coupling effects are suppressed beyond 2 T (inset Fig. 6), it is not surprising that quasi-2D excitations dominate the spectrum. We comment here that the CrX 3 materials have been considered as candidates for hosting topological magnons with Dirac cones in their (graphenetype) dispersion relations [22], but at finite applied fields a full gap is opened and the magnon dispersion has the quadratic form assumed in our model.
We stress again that, as shown in Fig. 1, these magnonic fluctuations are a significant source not of heat conduction but of phonon scattering. This situation is not generic, and indeed many of the systems reviewed in Sec. I provide examples of significant contributions to thermal conductivity from essentially 1D or 2D excitations in a 3D material [2][3][4][5][6][7][8][9]. In general one expects that, when such low-dimensional excitations are not mutually coherent between planes or chains with characteristic transverse coupling constant J ⊥ , they will act to destroy each others' coherence on a timescale governed by J ⊥ /J , where J is a characteristic energy scale (a gap or a stiffness) within the planes or chains. Thus, for a large system and non-negligible transverse coupling, the low-dimensional excitations would not contribute to thermal transport, as we observe for the 2D-coherent magnons of CrCl 3 , where J ⊥ /J ≡ J /3J 1/15, in region I. However, cases in which J is extremely small, the gap is very large, or other characteristics conspire to allow significant transport, have attracted attention precisely because they contradict this tendency towards incoherence.
The destructive effects of the 2D magnons on thermal conductivity are parameterized by the phonon scattering strength, λ(H, T ), in Eq. (4). To analyze λ(H, T ), in Fig. 7, we fit the thermal resistivity data, κ −1 (H ), for four of the five temperatures shown in Fig. 4. We find that an excellent description is obtained with a function λ(H, T ) ≡ λ(T ) completely independent of the field. Thus, in region I, the field dependence of κ (H, T ) is contained entirely within our straightforward assumptions about n mag (H, T ). Further, we observe that λ(T ) peaks in temperature close to T N , and thus we deduce that scattering by coherent but "only" 2D magnetic excitations has its strongest effect on the phonon thermal conductivity in the regime around T N itself; from the form of the magnetic specific heat [ Fig. 3(b)], it is no surprise that phonon scattering should be most efficient around this characteristic temperature. This leads to a suppression effect that includes the peak region (around T p ≈ 20 K) and becomes unusually large as the applied field is reduced below 2 T, where n mag rises strongly (Fig. 6). We note that, even in a van der Waals material, the phonons are in general much more 3D in nature than are the magnons once their interlayer correlations have been destroyed by the applied field, as a result of which good 2D magnons with no interlayer coherence are only damaging to phonons, and hence to thermal conduction.

B. Region II: dominant paramagnetic fluctuations
At higher temperatures, the magnetic fluctuations may no longer be regarded as well-defined 2D spin-wave excitations. As the energy scale of thermal fluctuations becomes large compared to that of the field, planar FM order is destroyed and the behavior of the spin system becomes paramagnetic. Despite this loss of coherent 2D magnons, magnetic scattering continues to play a significant role in controlling the thermal conductivity in region II, as Fig. 1 shows clearly. This effect must be a consequence of phonon scattering off randomly fluctuating free paramagnetic spins, whose relative density is also given by the fractional deviation of the magnetization from m s .
Although the dominant physics in region II is thermal fluctuations that flip individual spins against the magnetic field, the underlying spin interactions may by no means be neglected. We model the magnetization profile, m(T ), of CrCl 3 by a Weiss-field approach [41] in which the molecularfield term takes into account the FM correlations within the honeycomb planes. In this framework, where B S (y) is a Brillouin function of order S, S = 3/2, and g = 2. The solution for m(T ) is found by fixed-point iteration and, as shown in Fig. 8(a), an exceptionally good fit to the measured magnetization is obtained over the entire paramagnetic temperature range, T 40 K, for all fields below 7 T. The fitted constant B mol = 22.4 T 15 K is given once again by the FM in-plane coupling scale. In addition to providing smooth and easily computed curves at all measured field strengths, this model allows us to predict m(T ) with high confidence over the same temperature range at fields of μ 0 H = 9 and 18 T that are prohibitively high for SQUID magnetometer measurements (Fig. 8). The clear deviations between model and measurement at low fields and temperatures are consequences of the AF order, and are suppressed by fields in excess of 1 T (inset Fig. 6), to the point where deviations at the fields we can only model (μ 0 H = 9 and 18 T) are expected to be negligible.
To estimate the scattering strength, λ(H, T ) in Eq. (4), in region II, we observe first that when T → ∞ the magnetic degrees of freedom become entirely disordered, meaning that n mag → 1. Thus, to preserve the empirical observation that the field dependence of κ (T ) disappears at high T (Fig. 1), it is necessary that λ(H, T ) → 0 as T → ∞. On physical grounds, this must be the case because different phonon scattering mechanisms will overwhelm the magnetic contribution in Eq. (3), rendering τ −1 mag τ −1 0 . However, from the Weissfield estimate of m(H, T ) shown in Fig. 8, it is clear that this regime is reached only at our lowest measurement fields, and that n mag (H, T ) remains significantly less than unity over the available temperature range for all fields above 1 T.
From experiment, we have been able to obtain four highly accurate estimates of λ(T ) at temperatures up to T = 32 K, as shown in the inset of Fig. 7. Because the function λ(T ) is already well past its peak value, we use the logic of the previous paragraph to adopt test functions for the form of the vanishing of λ(T ) in region II, taking for specificity Gaussian, exponential, and power-law (Lorentzian) forms. Stated briefly, all forms of this functional tail give equally valid fits to the data throughout region II, and we pursue the quantitative aspects of this assertion in conjunction with the extraction of the phonon thermal conductivity, κ ph (T ), in the next subsection. The key qualitative point to be made here is that any H dependence of the functional tail is so weak that, once again, a completely adequate fit to all data is obtained by taking λ(H, T ) ≡ λ(T ), exactly as in region I.
Thus we have obtained the profound result that, to a very good approximation, the effect of a magnetic field on the thermal conductivity may be encoded entirely within the number of magnetic fluctuations, n mag (H, T ), acting to scatter the phonons transporting the heat, while their scattering efficiency is effectively H-independent.

C. Phonon thermal conductivity
In the preceding subsections we have shown that is it possible to obtain excellent descriptions of m(T ), and hence n mag (T ), in both regions I and II, as shown respectively in Figs. 6 and 8, by using minimal physical assumptions. Even without a quantitative matching of our two treatments, it is possible on this basis to obtain accurate fits of κ −1 (H, T ) at all fields and temperatures (Fig. 7) by using Eq. (4) with a field-independent scattering strength, λ(T ). The exercise of matching the models of Secs. IV A and IV B across the intermediate temperature regime would involve the extrapolation of either into a parameter range where it is explicitly no longer valid. However, the 2D spin-wave model applied in region I does account correctly for the temperature regime around the peak in κ (T ) for all fields in our measurement range, and thus this matching takes place only on the high side of the peak.
It is clear from Eq. (4) that κ (T ) is governed largely by the "continuous parameter" κ ph (T ), and as such that the most systematic matching of the two regimes would be ensured by obtaining an optimal estimate of this quantity. κ ph (T ), as the thermal conductivity due only to phonon contributions, is in principle obtained in region I as the H → ∞ limit, where the spin excitations have an infinite gap and n mag → 0, whence κ (T ) → κ ph (T ). In Sec. I, we summarized the difficulties in obtaining κ (T ) from a DC formalism [31,32], and in fact these are manifest even in obtaining κ ph (T ). Despite the many parameters, the form of the available relaxation-time terms remains highly constraining; in the present case, in the absence of a nonresonant spin-phonon suppression term, it is not possible to satisfy the low-and high-T limits simultaneously with an accurate estimate of the peak position, T p . An approach simpler than the full DC treatment is to proceed from the expression [1] where c ph (T ), the specific heat in the lattice sector, is known, v s is again an average phonon velocity, and a simple model can be constructed for the effective mean free path, eff , whose corresponding effective scattering time, τ eff (T ), may be connected to Eq. (3). However, once again it is not possible within such a simplified framework to obtain an accurate value of T p , and hence κ (H, T ) cannot be reproduced with any quantitative accuracy.
To ensure full generality of our treatment and to avoid the pitfalls inherent in adopting any approximate forms, we proceed instead directly from our data to compute the estimates κ est ph = κ (H, T )[1 + λ(T )n mag (H, T )] for each available field value. Our fits of κ −1 (H ) in region I, shown in Fig. 7, provide reliable estimates of κ ph at four discrete temperatures, which are shown as the green diamonds in Fig. 9. Inverting our κ (H, T ) data requires quantitative estimates of λ(T ) and a single function n mag (H, T ). For the latter we proceed, as shown in Fig. 8(b), by comparing our 2D spin-wave estimates at all low temperatures (solid lines) and our Weiss-field estimates at all high temperatures (dashed lines) with our m(T ) data at 1 and 5 T. As in Fig. 8(a), we also use our Weiss-field estimates at 9 and 18 T. By inspection, the upper temperature, T u (H ), beyond which the Weiss-field approach is quantitatively accurate, may be taken as T u (H ) = 40 K for all fields; in fact a constant T u (H ) would not be anticipated simply from the extent of region II in Fig. 5, and this result is actually a consequence of the additional effects of magnetic correlations above T N at low fields. The lower temperature, T l (H ), below which the 2D spin-wave result is quantitatively accurate, appears to be approximately 25, 30, and 35 K at 5, 9, and 18 T, values which track the edge of region I as represented in Fig. 5. At 1 T, the spin-wave approach is clearly no longer appropriate, because the tendency to AF order created by J is not fully suppressed ( Fig. 5; while we may use the m(T ) data here for quantitative purposes, we observe that the spin-wave result for a monolayer appears to be reliable to approximately 20 K. We note that these values of T l (H ) are not dissimilar to the values T max (H ) obtained from the peak in the magnetic specific heat (Sec. III C), and comment that this latter scale might indeed be anticipated as the upper limit to an ordered magnetic configuration on which to base a 2D spin-wave treatment. To bridge the shrinking gap between T l (H ) and T u (H ), we adopt the spline fits shown in Fig. 8(b). For λ(T ), as noted in Sec. IV B, we have tested three functions that reproduce the four data points spanning the peak in this quantity shown in the inset of Fig. 7; one has a Gaussian form, λ 1 (T ) = a 1 exp[−(T − T 1 ) 2 /2σ 2 ], one an exponential form, λ 2 (T ) = a 2 T b exp(−T /T 2 ), and one a Lorentzian form, The three functions differ only in the rate at which λ vanishes in the high-T regime and all three may be used to obtain consistent descriptions of the data with only minor differences in the value of the estimated κ ph (T ) at T > 40 K and in the temperature at which all datasets converge. Because of a mild but unexpected bulge around 40 K in our measured κ (T ) data (Fig. 1), which is most pronounced at 18 T and which we believe is not intrinsic to the phonon thermal conductivity, a rapid vanishing of λ(T ) cannot accommodate this feature. To avoid an ill-defined subtraction of this poorly understood contribution, it is most convenient to use the Lorentzian continuation of λ(T ), but we make no claim to have proven a physical underpinning for this form.
Inverting our κ (H, T ) data using these estimates of n mag (H, T ) and λ(T ) provides four different curves for κ est ph (T ) based on our measurements at 1, 5, 9, and 18 T. From Fig. 9, it is clear that all four estimates of κ ph (T ) are very similar across the full range of temperatures and that they converge completely in both the low-and high-T limits. In detail, the estimate based on our 1 T data, which are affected by the AF ordering tendencies, are something of an outlier. Otherwise, the three estimates based on 5, 9, and 18 T converge to high accuracy at all temperatures, with a maximum deviation of order 10% at T p . Because the minor bulge in our κ (T ) data around 40 K (Fig. 1) is strongest at 18 T, we judge the near-identical κ ph (T ) curves deduced from μ 0 H = 5 and 9 T to be the most representative and adopt these as our definitive result for the phonon thermal conductivity.
Fixing κ ph (T ) allows us unprecedented physical insight into the effects of phonon scattering by the spin fluctuations, and into the way in which these effects are suppressed by the magnetic field. In Fig. 10(a), we show the quantity λn mag , which is equivalent to the relative scattering rate τ 0 /τ mag . At 5, 9, and 18 T (dotted lines), λn mag displays a clear peak around 20 K that is suppressed from a value of nearly 1 by a factor of approximately 2 at each field step. For comparison, we show as solid lines the equivalent quantities for 0 and 1 T, which we have taken directly from our data because our methods for estimating n mag are not suitable at these low fields. We observe the same general form, but with the peak twice as strong again.
In Fig. 10(b), we show the suppression factor, f = (1 + λn mag ) −1 , corresponding to Fig. 10(a), which is an inverted function reaching a minimum of 1/3 at T p for a field of 0 T, where λn mag 2. Thus the quantitative conclusion concerning spin-fluctuation scattering of the phonon modes in CrCl 3 is that this mechanism is so effective around the Néel temperature of the AF ordered phase that 2/3 of the phonon contributions to the thermal conductivity are removed. The spin-scattering effect is not at all resonant, retaining a significant value over the entire range of temperatures from 0 to 100 K. Finally, in Fig. 10(c), we use our deduced values of κ ph (T ), n mag (H, T ), and λ(T ) to reconstruct our measured κ (H, T ) data for all fields and temperatures. Quantitatively excellent agreement is achieved in all cases, with deviations from 1% at 5 and 9 T to 15% around the peak at 1 T, where in any case our approximations are of limited validity. We stress that this is not a circular exercise, because all of the field dependence of κ (H, T ) is reproduced using only n mag (H, T ). We comment also that κ ph (T ) does lie significantly above our 18 T data; while this separation may be slightly exaggerated in the high-T regime [inset Fig. 10(c)] by our use of the Lorentzian continuation, the 20% difference around the peak position is a robust result that serves as a warning against assuming that high-field data must lie beyond the range of spin fluctuations.
It is clear that our models provide a quantitatively accurate fit to the thermal conductivity at all fields and temperatures with only one exception. This is the low-field plateau in the measured κ (T ) that arises due to the ordering transition, a feature that lies explicitly beyond our modeling. Thus we have shown that CrCl 3 provides an ideal system for modeling the magnetic fluctuations at finite fields and temperatures, and that the dramatic suppression of its thermal conductivity caused by these fluctuations can further be modelled by multiplying their density by a field-independent phonon scattering strength.

V. DISCUSSION AND CONCLUSION
We first reiterate what we have achieved in describing the thermal conductivity of a correlated magnetic insulator with significant spin-lattice coupling. We have shown that, beyond a very narrow regime of field and temperature hosting 3D magnetic order, the thermal conductivity is due entirely to phonons, but that the contributions of these phonons is subject to a suppression factor. This suppression is due entirely to scattering from spin fluctuations that are not coherent in 3D, and can be extremely strong (up to 65%), but can itself be suppressed to zero by a sufficiently strong magnetic field. This phenomenon can be captured quantitatively by expressing the suppression in terms of only two parameters, an independently determined number of active magnetic fluctuations and a dimensionless parameter for their scattering efficiency. This latter turns out to be entirely independent of the field, meaning that the scattering strength is dictated only by the temperature, while the full suppressive effects of the applied field are contained within the number of fluctuations.
Applying this very general framework, for CrCl 3 we model the fluctuation number, n mag (H, T ), by considering the two regimes of high and low H/T . In the former, magnetic fluctuations take the form of field-gapped 2D spin waves in the FM honeycomb planes and n mag is given by a conventional spinwave theory, but the complete lack of interlayer coherence means that these modes are not only incapable of transporting heat themselves but are destructive to the more 3D phonons that do so. In the latter, the system is paramagnetic and dominated by thermal fluctuations, but it is still field-polarized and thus the intrinsic intralayer FM interaction continues to play a role. The energy scale of this interaction, of order 15 K, is in fact fundamental to the thermodynamic and transport response of the system at all fields and temperatures, and its fingerprints are found in quantities ranging from the magnetic specific heat to the empirically determined phonon scatteringstrength parameter, λ(T ), which peaks in temperature around this value before falling to zero. With these ingredients, our formalism can reproduce, and indeed predict the form of, κ (T ) at different applied fields over the entire temperature range. Discrepancies between the data and the model appear only at low fields in the vicinity of T N , where n mag cannot describe adequately the population of fluctuating spins in the critical regime.
A key strength of our modeling procedure is its basis in the simple phonon-scattering form expressed in Eq. (4). Despite its phenomenological nature, our approach is both completely general and fully quantitative. It is not dependent on specific properties of the magnetic state of the system and, crucially, it is independent of the nature of the magnetic excitations, which in the quantum limit may not necessarily be conventional magnons, but topological ones or even only fractions of a single spin flip [3,33]. Our treatment is also insensitive to specific phonon scattering mechanisms, which as Sec. I makes clear can take many possible forms; in this context we note again that the spin-fluctuation scattering we consider is nonresonant and that these fluctuations are paramagnons (quasi-FM spin fluctuations in the paramagnetic regime) over the entire range of T and H. While we have kept our treatment independent of assumptions about the phonon thermal transport in order to make it fully quantitative (Sec. IV C), one may also ask whether the analysis could be improved, or otherwise brought into contact with conventional treatments based on the specific heat or the DC framework. While this is certainly possible in parts of the parameter space, we have not been able to find a global description for CrCl 3 within either approach, and on this basis would not expect this type of traditional formalism to be suitable for general magnetic materials. We stress again that, if the κ ph (T ) curve we extract from our data is regarded as a given, the number of "free" parameters in our modeling procedure can be argued to be zero, should one allow thatJ in region I and B mol in region II are given by the in-plane FM energy scale, 3J ≈ T N , and that the form of the high-T vanishing of λ(T ) is immaterial.
An essential aspect of our study is the issue of system dimensionality. While all magnetic insulators are 3D, in lowdimensional quantum magnets the regime of 3D behavior may be a very small part of the (H, T ) parameter space. In CrCl 3 , the stronger coupling in the honeycomb layers mandates a 2D treatment in the regime of large H/T , where the field destroys 3D correlations but does not damage the 2D FM correlations; by contrast, thermal fluctuations act to damage all correlations. The FM nature of the layers also has another unexpected consequence in that, although the CrX 3 systems are regarded structurally as van der Waals materials, featuring a very low cohesive energy for exfoliation [35], CrCl 3 does not show the conventional features of a 2D magnet. In the specific heat, the 3D ordering peak effectively coincides with the broad maximum characterizing the majority of the spinfluctuation processes, while the minimal anisotropy is also consistent with 3D magnetism. In more detail, the interlayer superexchange interaction, J ≈ 1 K, while not an insignificant fraction of the intralayer J = −5.25 K, is indeed rather smaller, and it is the FM nature of the in-plane order that allows J to "leverage" a T N scale of order 3|J|. We comment that layered magnetic materials are ubiquitous both in condensed matter and in the heterostructures being fabricated for spintronic functionalities, and hence our considerations can be expected to have far-reaching applicability.
While it is intuitively clear that the origin of the giant magnetoconductivity we observe lies in "strong spin-lattice coupling" [35], we stress that this is not merely another spin-phonon story. The qualitative difference in the present study is that we are measuring a transport property, meaning a property exclusively of the excitations in the spin and lattice sectors of the system. In this sense our focus is a specific and sensitive probe of a much less commonly studied aspect of magnetoelastic coupling, namely, the nature and scattering of these two sets of excitations over the complete (H, T ) parameter space.
One may nonetheless ask why, quantitatively, the scattering effect is so strong in CrCl 3 . Here we point to the possibility that the spin-phonon coupling can be relatively normal but J is in fact anomalously small. As noted in Sec. II, the FM J is a consequence of the near-90 • Cr-Cl-Cr bond angle enforced by edge-sharing CrCl 6 octahedra, and this type of interaction is far smaller in magnitude than comparable AF bonds at higher angles. In general, the interaction strength is very sensitive to the bond angle in this regime, implying that small phononic displacements may have strong relative effects; when normalized by the small J values, these magnetoelastic effects then appear in the conventional range. While first-principles calculations have been performed recently to accompany the experimental observation that the interlayer magnetic interaction strength, J , changes over a wide range in few-layer CrCl 3 samples [42], we are not aware of calculations investigating the lattice-sensitivity of the in-plane interaction, J, which could in principle be modulated by pressure in bulk samples or by substrate choice in few-layer heterostructured samples.
Returning to α-RuCl 3 , small Heisenberg interactions, J, are also a sought-after feature of candidate Kitaev materials. When these dominant magnetic effects are suppressed, the remaining terms in the anisotropic spin Hamiltonian are thought to give rise to fractional spin excitations (of Majorana [16,19,43] or generalized Majorana character [44]) and to spin-liquid ground states in the presence of an applied magnetic field [10,17,45]. It is clear [10] that the thermal conductivity measured in α-RuCl 3 shows strong spin-phonon scattering effects over the entire range of temperatures, but a detailed interpretation lies beyond the scope of a DC approach [21] and a theoretical analysis of phonon scattering by fractionalized spins is still awaited. We comment in passing that phonon coupling to chiral Majorana edge states has recently been invoked [46,47] as an ingredient essential for the interpretation of controversial thermal Hall conductivity data reported [48] for α-RuCl 3 at finite fields; however, we stress that the strong spin-phonon scattering observed in both α-RuCl 3 and CrCl 3 , and which we model here, involves the bulk spin excitations. We suggest that the more general, dataoriented approach we adopt for CrCl 3 may help to clarify the situation even in the absence of a microscopic discussion of phonon scattering by fractional spin excitations.
Returning to the higher chromium trihalides, CrBr 3 has been proposed [22] as a candidate for hosting topological magnons. CrI 3 is known [24] to present a situation where the bulk material has FM interlayer interactions, but the fewlayer form takes on a different interlayer structure and these interactions become AF. In a similar vein, it has been found very recently in CrCl 3 that the interlayer interaction remains AF [49] as the system thickness is reduced to two layers, while showing the dramatic increase noted above [42]. Efforts to include CrBr 3 in a systematic comparison are ongoing [50]. These results have drawn a great deal of attention with a view to fabricating highly controllable spintronic materials, possibly functioning with topologically protected information. While thermal conductivity measurements on few-or many-layer samples are not yet available, our results offer both a general framework for analyzing the different possible contributions to spin and thermal transport and a general warning concerning the need to take full account of spinphonon scattering effects.
In summary, we have investigated the thermal conductivity of the layered ferromagnet CrCl 3 over a wide range of temperature and magnetic field. We find a giant field-induced enhancement of the phononic contribution at all temperatures below 70 K, pointing towards a strong spin-fluctuation scattering effect. We construct an empirical model for the thermal conductivity by introducing a general framework based on two quantities that can be determined separately, the number of active spin-flip processes and their efficiency in scattering phonons. This formalism provides a quantitative description of our measured data at all fields and temperatures, has predictive power in unmeasured regions, and allows an accurate extraction of the purely phononic response. We anticipate that this approach will find wide application in interpreting the spin and thermal transport properties of many insulating magnetic materials, where spin-phonon scattering is a strong and unavoidable feature of the physics.