Introduction

Two-dimensional (2D) honeycomb lattices often exhibit topological properties, such as the quantum Hall effect, even in the absence of a magnetic field because of quantum confinement1. Graphene is such a system, hosting massless Dirac fermions at the Fermi surface. Bosonic systems can also have topologically nontrivial band structures that give rise to behaviors such as the spin Nernst effect2 and a thermal Hall effect3. Included in this list is CrX3 (X = Cl, Br, I), consisting of weakly bound van der Waals layers and considered to be magnetic analogues of graphene. Each layer is a 2D honeycomb lattice of magnetic Cr3+ ions with spin \(S = 3/2\). The ground state is either ferromagnetic (FM), as in CrBr3 and CrI3 with an out-of-plane spin orientation, or antiferromagnetic with an in-plane FM alignment that alternates in the perpendicular direction, as in the insulating CrCl3. All of the CrX3 compounds exhibit a structural phase transition driven by shifts in the layers from the high-temperature (HT) monoclinic \(C2/m\) phase to the low-temperature (LT) ABC rhombohedral \(R\bar 3\) phase (Fig. 1(a))4, though this transition is often broadened5. In CrCl3 the magnetic anisotropy is easy-plane while in CrBr3 and CrI3, there is an easy axis which is out-of-plane, and the interlayer coupling is stronger4,6. The interlayer magnetic coupling is stacking-dependent as demonstrated by a tenfold increase of AFM coupling strength reported in a thin \(C2/m\)-structure-stacked CrCl37. The edge-sharing octahedral coordination of Cr3+ with its 3d3 electronic configuration provides for superexchange spin-spin interactions. Unlike in the heavier halides, the spin-orbit coupling is expected to be much weaker in CrCl3.

Fig. 1: Structural details.
figure 1

a The crystal structure diagram of the \(R\bar 3\)-phase in CrCl3. b A honeycomb layer in CrCl3 with spins illustrated. c Powder diffraction data (intensity vs. lattice spacing d) collected on cooling from 300 to 5 K with Bragg peaks belonging to the two phases, labeled M-monoclinic and R-rhombohedral. The shaded region shows the diffuse scattering part of the intensity at 5 K. d Integrated intensity of the \(\left( {113} \right)_R\) and \(\left( {200} \right)_M/\left( {13\bar 1} \right)_M\) peaks plotted as a function of temperature, showing the progression of the transition on cooling from 300 K. Error bars smaller than symbols. e Integrated intensity of the magnetic \(\left( {00\frac{9}{2}} \right)_R\) peak, present below 15 K. Data collected on both cooling and warming. f Temperature-dependence of the a-axis lattice constant. A feature resembling negative thermal expansion is observed below ~50 K. g Interlayer spacing obtained from the position of the \(\left( {002} \right)_M/\left( {006} \right)_R\) Bragg peak. Data in (f) and (g) collected on cooling. Error bars represent one standard deviation.

Full size image

The magnetic transitions in CrCl3 are especially complex. From heat capacity, a single sharp peak corresponding to the FM transitions in CrBr38 and CrI39 is evident. For CrCl3, however, while a sharp peak is seen at the AFM transition at 14.1 K, an additional broad hump is present at 17.2 K5. This hump, as well as an inflection point observed in Faraday rotation data at 16.8 K10, has been interpreted as signaling the onset of a pseudo-FM phase having in-plane FM order but inter-plane disorder. Similar observations have been made for CrSBr, another layered compound with the same in-plane-FM and inter-plane-AFM order11. AC susceptibility measurements on CrCl3 provide additional evidence for two qualitatively different magnetic transitions; the real part of the magnetic susceptibility shows two peaks, located at 14.4 and 16.0 K, but the imaginary part only shows a single peak at 16.0 K12. Finally, we note that one study has reported the possibility of two or three transitions from low-field magnetic susceptibility measurements, as well as a hysteresis in the intensity of a magnetic peak measured by neutron diffraction on warming and cooling13.

Nuclear magnetic resonance (NMR)6,14,15,16,17 and inelastic tunneling spectroscopy18 have been used to probe the nature of the magnetic interactions in CrX3. In CrCl3, the nearest-neighbor in-plane exchange constant (J) is J = −0.90 meV from NMR in the range 0.4 ≤ T ≤ 8.1 K14, and J = −0.92 meV from inelastic tunneling spectroscopy at 0.3 K18. Interlayer interactions have also been probed by NMR, resulting in an out-of-plane AFM exchange constant of +0.003 meV (for a simplified model)14, and an AFM resonance19,20. Complementary to these techniques is inelastic neutron scattering (INS), as it can probe magnetic excitations across a wide range of wavevectors Q. CrBr321,22,23 and CrI324 have been probed by INS and both have been proposed to host topological magnons23,24. INS has also been used to investigate topological magnons in other honeycomb-lattice compounds, including CrSiTe325, CrGeTe325, and CoTiO326. Very recently, CrCl3 has been measured by INS27, but data from only a single temperature (4 K) has been reported. Studying a wide range of temperatures would yield insights on the nature of the magnetism in CrCl3.

We report elastic and INS measurements of the temperature dependence of the magnon and phonon excitations in CrCl3. Gapless Dirac magnons are observed at the lowest temperature where the dispersion intensity drops to zero at 4.5 meV, at the intersection of the optical and acoustic branches. On warming, the lower magnon branch shows a weak temperature dependence while the upper magnon branch shows a strong energy renormalization and lifetime broadening that is indicative of strong magnon-magnon interactions. Furthermore, an anomalous in-plane expansion of the honeycomb plane, along with signs of spin-phonon coupling, are observed below ~50 K; these are features of magnetoelastic effects in CrCl3 that suggest the onset of in-plane FM spin correlations well above the Néel temperature.

Results and discussion

Structure and lattice constant anomaly

Neutron scattering measurements were performed on a CrCl3 powder sample at the VISION instrument at Oak Ridge National Laboratory (ORNL), a time-of-flight spectrometer that simultaneously measures elastic and inelastic scattering data. The elastic intensity is shown in Fig. 1(c) as a function of temperature. Marked on the plot are Bragg peaks corresponding to the monoclinic \(C2/m\) and rhombohedral \(R\bar 3\) phases. It can be seen that the \(C2/m \to R\bar 3\) transition is broadened compared to that observed in single crystals5. Starting from 300 K, the structure is essentially single phase \(C2/m\), and all Bragg peaks can be indexed as such. Upon cooling, a gradual transition to the rhombohedral structure is observed, accompanied by diffuse scattering that most likely arises from a disordered sequence of HT- and LT-type layer stacking. Shown in Fig. 1(d) is the temperature dependence of the integrated intensity of the \(\left( {113} \right)_R\) peak and the overlapping \(\left( {200} \right)_M/\left( {13\bar 1} \right)_M\) peaks. Coexistence of both phases is observed over a wide temperature range in the powder sample starting at ~200 K and ending at ~25 K, that is distinctly different from the single crystal results reported earlier5. The transition proceeds gradually below 200 K until it stops around 25 K. Below ~15 K, the (0, 0, 9/2)R AFM Bragg peak is observed at d-spacing ~3.84 Å, consistent with the anticipated onset of AFM order at the Néel temperature (Fig. 1(e)). The lack of a discontinuity in the structural Bragg peak intensities near TN shows that the magnetic and structural phase transitions are not coupled.

An unusual behavior resembling negative thermal expansion (NTE) is observed in the a-axis lattice constant on cooling below ~50 K (Fig. 1(f)). We observed this behavior by tracking the position of the (300)R Bragg peak (see Supplementary Figure 2(a)). Data collected on powder and single-crystal CrCl3 on POWGEN (ORNL) and SPINS (NCNR), respectively, show the same temperature dependence in a; see Supplementary Fig. 3. In the out-of-plane direction, the interlayer spacing (obtained from the (002)M/(006)R peak position) plotted in Fig. 1(g) shows no clear anomaly, suggesting that the anomalous thermal expansion is confined in the honeycomb plane. The \(R\bar 3\)-phase unit-cell volume V can be calculated from a and the interlayer spacing, and we plot it in Supplementary Fig. 3(e); an anomalous flattening of V on cooling below ~50 K can be seen.

A similar in-plane NTE anomaly has been observed in other Cr3+ quasi-2D honeycomb magnets, namely CrBr328 and Cr2Ge2Te629, though in those materials the anomaly occurs near the magnetic ordering temperatures. (In CrI3, a small anomaly in the interlayer spacing has been reported9, but we are not aware of similar measurements of the a-axis lattice constant). Very recently, a theoretical study explained the NTE in CrBr3 as having magnons as a dominant contributor30. This study predicts similar behavior in CrCl3, though we note that their calculated relative change in the a-axis lattice constant between 5 and 15 K is roughly twice what we observe. The ilmenite compound NiTiO3, which shares the quasi-2D layered-honeycomb magnetism of CrCl3 but with 3D structural bonding, also exhibits a negative in-plane thermal expansion near its Néel temperature;31 however, in the structurally and magnetically similar compound CoTiO3, the in-plane thermal expansion is always positive, which was suggested to be due to a different electronic configuration for its Co2+ ions32.

We note that our observation of the in-plane NTE should not be affected by stacking changes. Since in-plane layer translations accompanying the LT- or HT-type stacking are nearly commensurate, Bragg peaks with HR and KR Miller indices divisible by 3 will be largely unaffected by variations in LT-type or HT-type stacking. (See Supplementary Note 1 for mathematical details). For illustration, Supplementary Fig. 1 shows the \(\left( {H0L} \right)\) plane of single-crystal X-ray scattering measurements; although the structural transition is broadened, the peaks along \(\left( {\bar 60L} \right)\) are undisturbed. Further cooling locks in the rhombohedral structure but with substantial diffuse scattering, evident underneath the Bragg peaks. Similar diffuse scattering has been noted for α-RuCl333.

Phonons

The in-plane NTE behavior coincides with changes in the phonon inelastic spectrum below ~50 K, most prominently near 26 meV (as first seen via Raman spectroscopy34). In our experiment, the INS intensity was measured along two narrow paths in (Q, ω) space, which are labeled the low-Q and high-Q trajectories in Fig. 3(a). These trajectories are determined by the fixed final neutron energy of 3.5 meV and the placement of the two sets of detectors at scattering angles of 45° and 135°. Intensity maps of \(\chi^{\prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) at temperatures from 5 to 300 K can be found in Supplementary Fig. 5. In Fig. 2(a), the dynamic susceptibility \(\chi^{\prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) along the high-Q trajectory is plotted at select temperatures. Several changes across the phonon spectra are observed with temperature, including the change at 26 meV indicated with an arrow. Shown in Fig. 2(b) is the fitted position of the 26 meV phonon peak as it rapidly shifts to higher energies upon cooling, starting around ~50 K and coinciding with the in-plane NTE. The 26 meV mode is a Cr-Cr stretching mode in the honeycomb plane (Fig. 2(d)).

Fig. 2: Phonons.
figure 2

a A plot of the dynamic susceptibility, \(\chi^{ \prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) as a function of temperature along the high-Q trajectory. Intensities have been smoothed; for clarity, error bars not shown. b Position of the 26 meV phonon feature from fitting. Error bars represent one standard deviation. c DFT-calculated phonon intensity along the high-Q trajectory for the \(R\bar 3\) structure with spins fixed to a certain ordering, including FM, in-plane AFM, or the inter-plane AFM observed in CrCl3. Inset depicts in-plane FM or AFM alignment. The features in the DFT intensity are about 7% lower in \(\hbar \omega\) than the corresponding features in (a). d Atomic movement of the 26 meV phonon mode.

Full size image

We performed density functional theory (DFT) calculations using three different magnetic structures: in-plane FM, in-plane AFM, or (the observed) inter-plane AFM, all of which have spins oriented in-plane. Overall, the data and calculations qualitatively agree, though the calculated spectra are lower overall by about 7%, and the 7 meV peak in the calculations is not present in the data. The spin-order dependence of the calculated spectra, though, does indicate the origin of the 26 meV phonon anomaly in the data. Slightly below 26 meV in the calculated spectra in Fig. 2(c), an anomalous decrease in the phonon energies is seen for in-plane AFM alignment relative to either in-plane FM or inter-plane AFM, showing that inter-plane order has little effect on the 26 meV phonon anomaly compared to the in-plane order. Since motion from the 26 meV mode involves Cr ion oscillations, this mode may induce frustration to the spin correlations. The differences in the relaxed lattice constants in the DFT calculations match the magnetoelastic coupling trends noted earlier, with a larger and c smaller in the FM alignment as compared to the in-plane AFM alignment, and little difference seen between FM and inter-plane AFM (see Supplementary Table 1). These calculations confirm that an increase in the in-plane spin-spin correlation on cooling, combined with the presence of spin-phonon coupling, is the cause of the phonon anomaly. It is remarkable that this substantial increase in spin correlations occurs well above the ordering temperature, on the order of 50 K ~3TN.

Magnetic excitations

The characteristic magnon spectrum of CrCl3 manifests as a drastic change in the low-energy \(\chi^{\prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) data between 5 and 50 K, as seen in Fig. 2(a). The temperature changes of the magnon intensity are shown in detail in Fig. 3(c–h), where a very strong temperature dependence is observed along both the low- and high-Q trajectories. A portion of the intensity is due to phonons, but judging from the 300 K data (Fig. 3(e, h)), the phonon intensity should be much weaker than the LT magnon intensity. (In detail, the magnetic \(\chi^{\prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) should decrease with temperature at high temperature. The phonon \(\chi^{\prime\prime}\left( {{{{\mathbf{Q}}}},\omega } \right)\), meanwhile, should (ideally) remain constant, which is consistent with the intensity “dip” near 10 meV being present at all temperatures and for both low- and high-Q. Thus, the 300 K data in Fig. 3(e, h) should present an upper bound on the phonon intensity, which should thus be negligible compared to the LT magnon intensity).

Fig. 3: Dispersion and temperature-dependence of magnetic excitations.
figure 3

a Powder-averaged \(S_ \bot \left( {{{{\mathbf{Q}}}},\omega } \right)\) simulated via SpinW, with low- and high-Q trajectories depicted in light blue. b Spin-wave dispersion in the \(\left( {HK0} \right)_R\) plane from fitted magnetic exchange coefficients. c–h Bose-factor-corrected inelastic neutron scattering intensity in the temperature ranges of (c, f) 5–20 K, (d, g) 20–100 K, and (c, h) 100–300 K, along the (c–e) low-Q and (f–h) high-Q trajectories. Insets of (e) and (h) show the uncorrected intensities of the data in those panels. All intensities smoothed, and error bars not included for clarity.

Full size image

To explain the observed magnon intensity, the spin waves were modeled using a Heisenberg Hamiltonian for a two-dimensional honeycomb lattice with FM nearest-neighbor (n.n.) interactions, which, analogous to graphene, features Dirac nodes at the K and K′ points of the Brillouin zone35, though we added coupling for 2nd- and 3rd-nearest in-plane neighbors as well. The dispersion is shown in Fig. 3(b). Due to the very weak interlayer magnetic coupling in CrCl314, we presume the dispersion is flat along the L-direction. (In Supplementary Fig. 4, we compare simulated magnon intensity for an interlayer coupling (per nearest interlayer bond) of +0.04 meV for the \(R\bar 3\) and \(C2/m\) structures to zero interlayer coupling, and find little difference. The comparison suggests that a relatively large \(J_L\) value of 0.04 meV is not enough to qualitatively change the intensity, and that the influence of stacking disorder on the magnon intensity is negligible). The calculated powder-averaged spin correlation function \(S_ \bot \left( {Q,W} \right)\) is shown in Fig. 3(a), obtained via the SpinW package in MATLAB36. The paths of the low- and high-Q trajectories are plotted as light blue curves, with the separation corresponding to the full-width-at-half-maximum of the energy resolution. The integrated intensity along the dispersion’s maximum, saddle points, and Dirac nodes leave sharply-defined features in the intensity along the low- and high-Q trajectories which have Q-independent energies (due to the flat dispersion along L.) From the energies of the Dirac node, saddle points, and dispersion maximum, we extracted magnetic coupling constants of \(J_1 = - 0.940\left( 5 \right)\) meV, \(J_2 = - 0.035\left( 2 \right)\) meV, and \(J_3 = + 0.070\left( 4 \right)\) meV. The uncertainties (from least-squares fitting) are likely underestimated, given the imperfect fit of the model to the data (see below); nevertheless, the obtained J1 values roughly agree with estimates from NMR14 and inelastic tunneling spectroscopy18.

The Dirac point manifests as a sharp cusp near 4.5 meV Fig. 4(a), where we show the 5 K inelastic intensity with the 300 K data subtracted to account for phonon intensity. Some magnetic intensity is also present at 300 K, especially at low-Q, but from Fig. 3(e, h) it appears to be be weak and featureless. The Dirac point is nestled between the upper and lower magnon dispersion branches. The red line in Fig. 4(a) shows the simulated intensity. Note that the Dirac node is gapless (a finding also reported in a recent INS study27), suggesting the absence of Kitaev or Dzyaloshinskii-Moriya (DM) interactions involving a DM vector component perpendicular to the c-axis. The dispersion saddle points result in peaks near 2.5 and 6 meV. The band maximum manifests as a steep drop-off in intensity near 8 meV. Although the data agree reasonably well with the J1-J2-J3 Heisenberg model, there are discrepancies, most notably a model peak near 7 meV at low Q which is not seen in our data. This peak corresponds to the crossing of the low-Q trajectory with a boundary (see Fig. 3(a)) corresponding to where new magnon branches become accessible as Q changes at constant \(\hbar \omega\). We speculate that, even at 5 K, magnon-magnon interactions might cause disorder in the spin orientations that suppress this peak. More studies are needed to explore this issue further.

Fig. 4: Analysis of spin waves.
figure 4

a Comparison of 5 K inelastic data with calculated \(\chi^{\prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) along low- and high-Q trajectories. The 5 K data has the 300 K data subtracted to eliminate the contribution from the 7 to 10 meV phonon feature, which can be seen in, e.g., the high-Q data in Fig. 3(h). Data are unsmoothed. b Intensities near 1 meV, the Dirac node energy, and the middle of the upper branch (“max pt”) for the low-Q trajectory. The Dirac node energies and the position and energy of the midpoint of the upper branch hump were determined by fits to a function of connected line segments. The intensities near the Dirac node and 1 meV were taken as an average within ±0.1 meV of these locations. c The center-of-mass (c.o.m.) position of the upper branch intensity plotted as a function of temperature. The black curve shows the expected first-order renormalization \(\hbar \omega = \hbar \omega _0\left( {1 - \alpha _1T^2} \right)\) of the honeycomb ferromagnet with nearest-neighbor interactions, where \(\alpha _1 = \pi /\left( {24\sqrt 3 J^2S^3} \right)\)35. Error bars represent one standard deviation.

Full size image

On warming, the dispersion features soften and broaden (as seen in Fig. 3(c–h)), but many features remain well-defined above TN. The “corner” near 6 meV at low-Q remains visible at 20 K, and a “hump” of upper branch intensity is still present at 60 K. The valley near 4.5 meV shifts to lower frequency due to magnon energy renormalization, indicating that spectral weight above and below the Dirac node fills in the surrounding energy range. The upper edge near 8 meV, already broader than resolution at 5 K (Fig. 4(a)), becomes markedly less steep on warming. Notably, the peak near 2.5 meV, especially large at high-Q, vanishes around TN after shrinking steadily with temperature.

To characterize the temperature changes in the magnon intensity, the integrated intensity at three points of the 5 K data in Fig. 4(a) (marked with arrows) are plotted as a function of temperature in Fig. 4(b). These points are at 1 meV, the Dirac node, and “max pt” (the midpoint between the two corners of the upper branch intensity, which are ~6 and 7.5 meV in Fig. 4(a)), with the latter two points obtained via a fit to a function of connected line segments (see Supplementary Note 7 and Supplementary Fig. 6). The data were integrated within ±0.1 meV of these \(\hbar \omega\) values. Starting with the Dirac node, the intensity is ~zero at 4.5 meV, but gradually increases because of intensity spilling over from the upper branch due to the magnon energy renormalization. In contrast, the intensities at the other two points, representative of the lower and upper magnon branches, decrease with temperature, though, below TN, the lower branch changes more gradually than the upper branch. In Fig. 4(c), we show the center-of-mass (c.o.m.) position of the upper branch intensity. It turns out that the higher energy branch softens more rapidly than can be accounted for by the first-order magnon-magnon renormalization for a honeycomb ferromagnet35, shown as the black curve. Thus, the magnons in CrCl3 exhibit significant broadening in energy, and a shift to lower energies, as shown by the ~2 meV (~25 %) decrease in the c.o.m. position from 5 to 50 K. As noted, even at 5 K we observe a broadening that is greater than resolution, in contrast to a recent study (on an instrument with lesser resolution) which was not able to observe such broadening at 4 K27. Our data at a multitude of additional temperatures further emphasize the substantial broadening of the magnon features on warming.

At high temperature, we expect the spin waves to become overdamped, leading to quasi-elastic scattering with a Lorentzian lineshape centered at 0 meV, and no change with temperature assuming neither the magnon decay rate or Q-coefficients change37. In the insets of Fig. 3(e, h), we see that the intensity (without Bose factor correction) is roughly constant from 100 to 300 K, especially for the low-Q trajectory where phonons have less contribution (though at high-Q, some phonon-induced intensity increase with temperature can be seen). This temperature-independent intensity is consistent with expectations. In Supplementary Fig. 7, we show that the overall intensity (\(\chi^{\prime\prime} \left( {{{{\mathbf{Q}}}},\omega } \right)\) integrated from 1 to 10 meV) decreases more rapidly at high-Q than low-Q, especially at low temperature.

The onset of spin waves occurs on the same temperature scale as the in-plane FM spin correlations, as indicated via the anomalies in the a-axis lattice constant and the phonon feature; these behaviors all occur gradually on cooling below ~50 K, 10s of K higher than TN. In contrast, CrBr3 exhibits abrupt changes at the Curie temperature in the spin-phonon coupling and the lattice constants28, and these differences may be due to differing magnetic anisotropy. Anisotropy can play an important role in the ordering behavior of 2D systems; the 2D isotropic Heisenberg model has no long-range order above \(T \,> \,0\) for short-range interactions38, in contrast to the transition present in the 2D Ising model39. The different spin orientations of the chromium trihalides hints at a different magnetic anisotropy for CrCl3 than for CrBr3 and CrI3. Though estimates from inelastic tunneling spectroscopy suggest that the anisotropy is on the order of a few percent or less18, even a slight asymmetry may influence magnetic behavior, such as the presence of a Kosterlitz–Thouless transition in K2CuF440. Studying the magnetic excitations of CrCl3 in detail may elucidate the role of anisotropy in the magnetism of quasi-2D materials.

A diverse array of materials have been investigated for their topological properties; CrBr323, CrI324, and CrXTe3 (X = Si, Ge)25 have been reported to host topological magnons In contrast, our CrCl3 data show a sharp cusp at the Dirac point in powder-averaged INS intensity, indicating no gap within the precision of our experiment. CrCl3 thus provides a foil to the other chromium trihalides, and may facilitate understanding of spin dynamics across the CrX3 system. Meanwhile, the remarkable temperature dependence of CrCl3, in which the magnon intensity has only subtle changes across TN and continues to have identifiable upper branch intensity tens of Kelvin above TN, suggests that the temperature dependence of CrBr3 and CrI3 should be investigated in more detail, especially above \(T_C\). For CrI3, on approaching \(T_C\), the magnetic anisotropy gap appears to vanish while the acoustic magnon stiffness remains substantial24, but data well above \(T_C\) has not been reported and would be informative. Further study should clarify the extent to which the behavior of CrCl3 is unusual.

In conclusion, from elastic and INS measurements, we have observed changes in the structural, vibrational, and magnon intensity as a function of temperature that share a similar temperature scale. At 5 K, the magnon intensity shows a sharp cusp indicating a dispersion with gapless Dirac magnons. On warming, a substantial energy renormalization is seen that occurs gradually on warming, with little discontinuity at TN and upper branch intensity still visible tens of Kelvin above TN. Below ~50 K, a planar NTE and an increase in the energy of a phonon feature near 26 meV can be seen. These anomalies are explained via DFT calculations as resulting from spin-phonon and magnetoelastic coupling combined with the increase in in-plane spin correlations on cooling. Altogether, our results suggest a persistence of in-plane spin correlations to the order of ~3 TN, as seen from the behavior of the a-axis lattice constant, phonons, and magnons.

Methods

Neutron scattering

Most data were collected on VISION41, an indirect-geometry time-of-flight neutron spectrometer at the Spallation Neutron Source (SNS) at ORNL. The final neutron energy was fixed at 3.5 meV. Inelastic data were taken on two detector banks at scattering angles of 45 and 135°, while six elastic scattering detector banks were located at 90°. The CrCl3 powder was poured into a vanadium can. Background from an empty-can scan was subtracted from all data. Measurements were also taken at the time-of-flight instrument POWGEN42 (at the SNS at ORNL) and on the triple-axis spectrometer SPINS (at the NIST Center for Neutron Research).

Sample details

The powder samples measured in VISION and POWGEN consisted of as-purchased CrCl3 powder from Alfa Aesar (99.9% purity, metals basis), consisting of mm-sized flakes. The single crystal measured on SPINS was grown via the vapor transport growth procedure detailed in ref. 5.

Density functional theory calculations

Spin-polarized Density Functional Theory (DFT) calculations of CrCl3 (in its \(C2/m\) and \(R\bar 3\) phases and with FM in-plane AFM, and inter-plane AFM configurations) were performed using the Vienna Ab initio Simulation Package43. The calculation used the Projector Augmented Wave method44,45 to describe the effects of core electrons, and the Perdew–Burke–Ernzerhof (PBE)46 implementation of the Generalized Gradient Approximation (GGA) for the exchange-correlation functional. The energy cutoff was 600 eV for the plane-wave basis of the valence electrons. The lattice parameters and atomic coordinates from ref. 47 were used as the initial structure, and they were then fully relaxed to minimize the potential energy. The electronic structure was calculated on a Γ-centered mesh (9 × 5 × 9 for \(C2/m\) and 9 × 9 × 3 for \(R\bar 3\)). The total energy tolerance for electronic energy minimization was 10−8 eV, and for structure optimization it was 10−7 eV. The maximum interatomic force after relaxation was below 0.001 eV Å−1. The optB86b-vdW functional48 for dispersion corrections was applied, and a Hubbard U term of 3.7 eV49 was applied to account for the localized 3d orbitals of Cr. A supercell (2 × 1 × 2 for \(C2/m\) and 2 × 2 × 1 for \(R\bar 3\)) was created for phonon calculations. The interatomic force constants were calculated by Density Functional Perturbation Theory, and the vibrational eigenfrequencies and modes were then calculated using phonopy50. The OCLIMAX software51 was used to convert the DFT-calculated phonon results to the simulated INS spectra.

A standard deviation of statistical uncertainty is denoted with parentheses in the last digit(s) for numbers, and with error bars for plots.