Field evolution of magnetic phases and spin dynamics in the honeycomb lattice magnet Na₂Co₂TeO₆: ²³Na NMR study

Figure 2 shows the temperature and field dependences of the magnetic susceptibility $M/B$ under the ZFC condition. No essential difference was observed between the ZFC and FC conditions except a minor difference at 1 T below $T_{N}$ (not shown) which may be attributed to compensation behavior of ferrimagnetism [45]. The results are in good agreement with those reported on polycrystalline samples [36, 37, 38] but are slightly different in the absolute magnitude from that of a single crystal because of powder averaging of the anisotropic susceptibility [45, 46]. The Néel temperature and its field variation also agree with those in the previous reports. The susceptibility is essentially field independent in the PM phase, contrasting to the behavior in $\alpha$-RuCl₃ where it is enhanced by field at temperatures $T\lesssim 3T_{N}$ [60].

Figure 3(a) shows the temperature ($T$) variation of the ²³Na NMR spectrum taken at a frequency $\nu_{0}=33.790$ MHz. The signal intensity is plotted against the field shift $B_{0}-B$ relative to the reference field $B_{0}=\nu_{0}/\gamma$ where $B$ is the external field. The spectra at high $T$ in the PM phase exhibit a quadrupolar-split powder pattern characterized by asymmetric electric field gradients (EFGs) at the Na site [61]. This is in clear contrast to the previous reports on single crystals, none of which observed quadrupolar satellites [44, 49]. The $T$-dependent magnetic broadening is also obvious (see below). From singular positions of the satellites, the principal values of the EFG tensor are determined as ${}^{23}\nu_{\mathrm{Q}}=1.68$ MHz and $\eta=0.49$ which are almost $T$ independent. We observed neither line splitting due to the presence of several crystallographic Na sites nor severe broadening or smearing of quardupolar satellites due to disorder [36, 38, 37]. This suggests preferential occupation of some Na site and/or a tendency of atomic ordering rather than strong disorder in Na layers.

We measured the NMR spectrum at various frequencies and temperatures and determined the temperature and field dependences of the hyperfine field at the Na sites and the line width of the spectrum. In the following, we label the spectra and the related quantities by the reference field $B_{0}$. Figure 4(a) shows the $T$ dependences of the line width (fwhm) $\Delta B$ taken at various $B_{0}$’s. $\Delta B$ is also plotted in Fig. 4(b) as a function of uniform magnetization $M$ per Co atom with $T$ the implicit parameter. It was found that $\Delta B$ in the PM phase above 50 K is proportional to $M$ in a wide range of $B_{0}$. The proportionality constant $\Delta B/M=108\mathrm{~{}mT}/\mu_{B}$ is comparable to the root mean square of the principal values of the dipolar field tensor of $81\mathrm{~{}mT}/\mu_{B}$ calculated for the most occupied ($\sim$70%) $12i$ site for Na. The line width in the PM phase is thus dominated by anisotropic dipolar coupling between Na and Co.

As shown in Figs. 3(a) and 4(a), sudden increase of the line width was observed on entering the AFM phase at relatively low fields $B_{0}\leq 5$ T. The low-field line widths exhibit order-parameter-like $T$ dependence, saturating at low $T$ with a slightly field-dependent value. The spectrum is broadened almost symmetrically about the reference field $B=B_{0}$, and $\Delta B$ no longer scales with $M$ [Fig. 4(b)]. These observations demonstrate the appearance of a staggered hyperfine field at the Na site of which magnitude is much smaller than the external field. The anomaly in $\Delta B$ at $T_{N}$ is obscured at high fields $B_{0}\geq 7$ T owing to magnetic broadening present already in the PM phase. The high-field line widths saturate in the low-$T$ limit to a value somewhat larger than the low-field value, possibly related to a reorientation of ordered moments with the external field.

It should be emphasized that at low $T$ the scaling between $\Delta B$ and $M$ breaks down even at the highest field $B_{0}=8.9$ T. This is also true for the hyperfine field $B_{\mathrm{hf,pk}}$ at the peak position of the spectrum as will be shown below [Fig. 4(d)]. If there appeared a field-induced spin-disordered phase as suggested in Ref. 46, both $\Delta B$ and $B_{\mathrm{hf,pk}}$ should scale linearly with $M$ down to the lowest temperature because of the disappearance of a staggered component of the hyperfine field. The upward deviation of $\Delta B$ from the $\Delta B$-$M$ scaling means that the spectrum is much broader than is expected from the uniform magnetization, demonstrating a significant contribution of the staggered hyperfine field to $\Delta B$. Our results thus point to the absence of a field-induced disordered phase up to $B\sim 9$ T consistent with a report that the in-plane critical field is around 10 T [48].

Figure 3(b) shows the field evolution of the NMR spectrum at 4.2 K. The line shape as well as the line width does not change much up to $B_{0}=8.9$ T, which confirms the persistence of magnetic LRO. A shift of the peak due to slight asymmetry in the line shape is discernible at $B_{0}\geq 2$ T. The peak shifts to a positive side of $B_{0}-B$ above $B_{0}=7$ T, which may be associated with a magnetization jump around 6 T for the field applied perpendicular to the zigzag chains [39, 45].

The magnetic shift $K$ is an important measure of the local static spin susceptibility of a magnetic ion. For nuclei with $I\geq 1$, care must be taken to correct a contribution of the quadrupole interaction to the shift of the NMR line [61]. We measured the field dependence of the relative line shift ${(B_{0}-B)/B}$ and confirmed that the second-order quadrupolar shift is negligible compared with the magnetic shift for $B_{0}\geq 3$ T. The field shift $B_{0}-B$ near the peak of the spectrum can then be taken as the hyperfine field $B_{\mathrm{hf}}$ at the Na sites from which the magnetic shift is defined as $K=B_{\mathrm{hf}}/B$.

Shown in Fig. 4(c) is the $T$ dependence of $B_{\mathrm{hf,pk}}=B_{0}-B_{\mathrm{pk}}$ for various $B_{0}$’s where $B_{\mathrm{pk}}$ is the external field at the peak position of the spectrum. Figure 4(d) shows a plot of $B_{\mathrm{hf,pk}}$ against $M$ alternative to the conventional $K$-$\chi$ plot ¹¹1This plot is useful in determining the hyperfine coupling constant when magnetization shows non-linear increase at high fields so that the susceptibility $\chi$ defined by $M/B$ depends on a magnetic field.. The linear relation $B_{\mathrm{hf,pk}}=a_{\mathrm{hf}}M$ holds in the PM phase in a wide range of $B_{0}$. The slope $a_{\mathrm{hf}}=23.7\mathrm{~{}mT}/\mu_{B}$ is by definition the hyperfine coupling constant between Na and Co. The small and positive value of $a_{\mathrm{hf}}$ suggests that the Na-$3s$ orbital is polarized a little due to spin transfer from the neighboring Co atoms. Breakdown of the scaling between $B_{\mathrm{hf,pk}}$ and $M$ in the AFM phase is signaled by deflection of the plot from the straight line. The deviation of the low-$T$ data at the highest field $B_{0}=8.9$ T from the scaling confirms the absence of a spin-disordered phase up to $B\sim 9$ T.

Another interesting behavior in the AFM phase is a field-dependent shift of $B_{\mathrm{hf,pk}}$ accompanying the sign change [Fig. 4(c)]. As noted above, this is due to asymmetry in the line shape that reflects the distribution of staggered hyperfine fields. The negative shift of the peak at low fields may be associated with domains of net magnetic moment directed opposite to the external field which exist in ferrimagnets showing compensation behavior. The positive shift at high fields is probably caused by a decrease in the number of such domains as well as an increasing contribution of uniform magnetization to the hyperfine field.

Integrated intensity of the NMR spectrum is one of the essential quantities to be measured in frustrated magnets because it is highly sensitive to slow dynamics of which existence is manifested by the loss of signal intensity, or wipeout [63, 64, 65]. Recently, the slow dynamics of Co spins in the AFM phase has been suggested from partial wipeout of the ²³Na NMR signal below $T_{N}$ [44]. To avoid an artifact of finite separation $\tau$ between rf exciting pulses in determining the integrated intensity, we measured spin-echo decay at the peak position of the spectrum and corrected the intensity by extrapolating it to $\tau=0$. The results are shown in Fig. 5 where we plotted the corrected intensity multiplied by $T$, a factor arising from nuclear paramagnetism, as a function of $T$. Unlike the previous report, we observed no anomaly around and below $T_{N}$; the integrated intensity is essentially $T$ independent down to 4.2 K. This clearly shows the absence of slow dynamics in the measured $T$ range ²²2We detected a dip in the uncorrected intensity (typically taken at $\tau=40\mathrm{~{}\mu s}$) around $T_{N}$ as reported in Ref. 44 which is recovered by the extrapolation process. The apparent decrease in intensity around $T_{N}$ reported in Ref. 44 could be caused by an insufficient correction of shortening of spin-spin relaxation time $T_{2}$ accompanied by critical slowing down. The unrecovered signal loss below $T_{N}$ might also be an artifact of insufficient integration range resulting from the large spectral broadening due to magnetic LRO.

Figure 6(d) shows examples of the recovery of ²³Na spin-echo intensity $M(t)$. The data are well fitted by Eq. (1), which demonstrates reliability of the analyses. The $T$ and $B$ dependences of $1/T_{1}$ at the Na sites are shown in Fig. 6(a), and the low-$T$ close-up in the form of $1/T_{1}T$ in Fig. 6(b). The ²³Na $1/T_{1}$ in Na₂Ni₂TeO₆ with a similar Néel temperature $T_{N}\approx 26$ K but with a slightly different stacking sequence of honeycomb layers are shown for comparison [67]. The ²³Na $1/T_{1}$ in Na₂Co₂TeO₆ is characterized by a strong variation not only with temperature but with magnetic field. This makes a striking contrast to the macroscopic susceptibility which depends hardly on magnetic field in the PM phase up to 9 T [Fig. 2]. The four characteristic temperature regions are identified: Region I, above 50$-$60 K ($\sim 2T_{N}$) where $1/T_{1}$ is nearly $B$ independent and the $T$ dependence is relatively weak; Region II, $T_{N}<T\lesssim 2T_{N}$ where $1/T_{1}$ increases toward $T_{N}$ and the $B$ dependence becomes noticeable; Region III, $T_{N}/2\lesssim T<T_{N}$ where $1/T_{1}$ decreases with decreasing $T$ but remains enhanced, still exhibiting strong field evolution; and Region IV, $T\lesssim T_{N}/2$ where a rapid decrease of $1/T_{1}$ is observed regardless of the field. On the $B$ dependence, $1/T_{1}$ shows a contrasting response to magnetic fields above and below $T_{N}$; $1/T_{1}$ is suppressed by field above $T_{N}$ whereas it is enhanced below $T_{N}$.

Let us inspect first the $T$ dependence of $1/T_{1}$ at the lowest field. $1/T_{1}$ measured at $B=1$ T depends only weakly on $T$ above about 60 K, approaching an almost $T$- and $B$-independent value of $1/T_{1\infty}\sim 100\mathrm{~{}s^{-1}}$ as usually observed in magnets with exchange-coupled local moments. On decreasing $T$ across 60 K, $1/T_{1}$ starts to increase due to the development of short-range spin correlations consistent with the neutron diffuse scattering [38]. This is followed by a divergent increase of $1/T_{1}$ on approaching $T_{N}$ below about 30 K, an indication of three-dimensional (3D) critical slowing down. $1/T_{1}$ takes a maximum at $T_{N}\approx 26.5$ K and then decreases rapidly on cooling. However, the decrease of $1/T_{1}$ is rather gradual compared with that in a conventional antiferromagnet in which $1/T_{1}$ often decreases many orders of magnitude not far below $T_{N}$ as observed in Na₂Ni₂TeO₆. $1/T_{1}$ remains in the same order as in the PM phase, which indicates residual low-energy spin excitations in Region III. This is closely related to the result of single-crystal INS that the low-energy spectral weight survives down to $\sim$14 K without forming a spin-wave mode and a gap in the excitation spectrum [44]. The absence of an upturn of thermal conductivity on entering the AFM phase [48] may have a close connection to the residual low-energy spin excitations. On further cooling below 13$-$15 K to Region IV, $1/T_{1}$ decreases rapidly over two decades. The $T$ dependence of $1/T_{1}$ is approximated by a power law $1/T_{1}\propto T^{n}$ with $n\approx 5$ predicted for the three-magnon process at $T\gg\Delta$ rather than an activation law expected at $T\ll\Delta$ where $\Delta$ is the gap in the spin-wave spectrum [68].

Our result at $B=1$ T is in good agreement with the previous report of $1/T_{1}$ for a single crystal by Chen et al. taken at $B=0.75$ T with $\mathbf{B}\parallel a^{*}$ from 15 to 60 K [44]. Minor differences in the absolute magnitude of $1/T_{1}$ may be attributed to our use of the stretched exponential function Eq. (1) as well as of polycrystals ³³3The stretched exponential function used in Ref. 44 is slightly different from conventional ones: the exponential terms are like $\exp[-6(t/T_{1})^{\beta}]$ rather than $\exp[-(6t/T_{1})^{\beta}]$. If this form of stretched exponentials is applied to our recovery, the resulting $1/T_{1}$ agrees quantitatively with their $1/T_{1}$.. A more extensive single-crystal study has been reported by Lee et al. who measured both $1/T_{1}$ and $1/T_{2}$ at $B\sim 3.1$ T with $\mathbf{B}\parallel c$ and $\perp c$ up to room temperature [49]. The $T$ dependence of their $1/T_{1}$ above $\sim$10 K is qualitatively similar to ours at 3 T, but the absolute values are larger by a factor of $\sim$2. The discrepancy becomes progressively greater below $\sim$10 K in Region IV, leading to a moderate $T$ variation of their $1/T_{1}$ described by a power law with a smaller exponent $n\sim 3$.

In the AFM phase where the powder NMR spectrum is largely broadened due to the staggered hyperfine field, we measured $1/T_{1}$ at several positions of the spectrum other than the peak position and found that the $1/T_{1}$ differs by at most 10% in parallel with the result of $1/T_{1}$ for the split peaks in Ref. 49. We also performed the inverse Laplace transform analysis of magnetization recovery [70, 71] in order to see whether the stretched exponential analysis captures the true distribution of $1/T_{1}$. The distribution function of $1/T_{1}$ similar to that expected for the stretched exponential recovery [72, 73] was obtained as detailed in the Appendix. This justifies our phenomenological analysis of $1/T_{1}$ using Eq. (1) and rules out a possibility of the inequivalent Na sites showing a distinct $T$ variation of $1/T_{1}$ from the one displayed in Fig. 6.

The origin of the large discrepancy between our $1/T_{1}$ and that reported in Ref. 49 especially in the behavior in Region IV is unclear. This is partly because the authors of Ref. 49 did not give fundamental information on the recovery such as the $T$ dependence of $\beta$ to be compared with ours and the recovery itself from which $1/T_{1}$ is extracted. A possible origin of the discrepancy is a different degree of atomic disorder in Na layers between the samples; their single-crystal sample may have stronger disorder than our polycrystals because quadrupolar satellites are missing in their NMR spectrum ⁴⁴4Vanishingly small quadrupolar splitting or negligible quadrupole coupling seems improbable for their missing of quadrupolar satellites because they can fit the magnetization recovery using the form for a quadrupolar-split center line.. Disorder in Na layers might change low-energy spin excitations by modulating interlayer coupling and/or by introducing bond randomness in Co layers, affecting $1/T_{1}$ effectively at low $T$ where the intrinsic $1/T_{1}$ falls off due to magnetic LRO. We should also point out that in the PM phase the Redfield contribution ($T_{1}$ process) to $1/T_{2}$ for $\mathbf{B}\parallel c$ evaluated from their $1/T_{1}$ for $\mathbf{B}\parallel c$ and $\perp c$ exceeds the observed $1/T_{2}$ for $\mathbf{B}\parallel c$, which is impossible for the quadrupolar-split center line [75, 76]. To resolve the conflicts, single-crystal NMR measurements are under way and will be reported in a future publication.

Let us go back to our results at higher fields. Increasing magnetic field suppresses $1/T_{1}$ above $T_{N}$, most drastically around $T_{N}$ where the peak gets broadened and is shifted to lower $T$. This suggests field-induced suppression of short-range spin correlations. In contrast, $1/T_{1}$ is enhanced by field below $T_{N}$ to exhibit a shoulder-like anomaly around 10 K. The anomaly appears as a broad hump in a plot of $1/T_{1}T$ and is most pronounced at 6$-$7 T [Fig. 6(b)], which indicates field enhancement of the spectral weight remaining at low energies. The field enhancement of $1/T_{1}$ cannot be ascribed to defect spin fluctuations because in that case $1/T_{1}$ would be suppressed by magnetic fields [53, 77]. On the other hand, the magnetic field does not affect much the $T$ dependence of $1/T_{1}$ in Region IV; $1/T_{1}$ shows an approximate $T^{5}$ law even at 9 T. There is no indication of strong field enhancement of the low-energy excitations. This corroborates the absence of a spin-disordered phase up to 9 T.

The rounding of the peak of $1/T_{1}$ around $T_{N}$ and the appearance of a broad hump in $1/T_{1}T$ around 10 K in magnetic fields resemble the behavior of magnetic specific heat $C_{m}/T$ [45, 46]. This suggests that $1/T_{1}$ probes excitations governing the magnetic specific heat. The resemblance between $1/T_{1}T$ and $C_{m}/T$ also suggests that the rounded peak of $1/T_{1}$ in magnetic fields is not an artifact due to the use of powder samples but an intrinsic property of this compound ⁵⁵5Since the magnetic response of Na₂Co₂TeO₆ is insensitive to the field applied perpendicular to the honeycomb planes, it seems reasonable to consider that the field dependence measured in powders is governed by in-plane field responses. This is further supported by the fact that in powders a grain with the field lying in the honeycomb plane is found more frequently than a grain with the field normal to the plane, because the probability of a field making an angle $\theta$ with respect to the direction normal to the plane is proportional to $\sin\theta$..

The stretching exponent $\beta$ in Eq. (1) also includes valuable information on the spin dynamics. As shown in Fig. 6(c), $\beta$ is nearly $T$- and $B$-independent from 30 to 200 K and takes a value close to unity, which indicates a nearly uniform relaxation process. A rapid decrease of $\beta$ above 200 K implies the appearance of additional relaxation channels possibly related to the spin-orbit excited state lying 22$-$23 meV above the ground state [41, 42]. This provides an indirect support for the spin-orbital entangled state of Co²⁺ in Na₂Co₂TeO₆. $\beta$ decreases steeply below 30 K and takes a local minimum around $T_{N}$. A large distribution of $1/T_{1}$ close to $T_{N}$ may result from strong temperature and field-orientation dependence of $1/T_{1}$. While modest inhomogeneity of $1/T_{1}$ ($\beta\gtrsim 0.9$) is found above $\sim$10 K, inhomogeneous relaxation prevails very rapidly below $\sim$10 K. This suggests a qualitative change of the spin dynamics on entering Region IV.

The magnetic phase diagram is constructed from the $T$- and $B$-dependences of various quantities measured by NMR. The result is summarized in Fig. 7 together with a contour plot of $1/T_{1}T$ for comparison. $T_{N}$ may be determined in several ways: the temperature below which the $B_{\mathrm{hf,pk}}$-$M$ scaling breaks down [Fig. 4(d)], the temperature at which $1/T_{1}T$ takes a maximum, and the temperature at which $\beta$ takes a local minimum. All of these agree within experimental accuracies. Determination of the boundary between Regions III and IV seems more difficult because of a gradual nature of the transition. We tentatively adopt the temperature at which the derivative of $1/T_{1}T$ takes a local maximum. The field above which $B_{\mathrm{hf,pk}}$ shifts to a positive side [Fig. 3(b)] constitutes another boundary dividing low- and high-field AFM phases, although the small shift of $B_{\mathrm{hf,pk}}$ relative to the broad powder spectrum leads to large errors in the boundary field. The phase diagram is in good agreement with that determined from the macroscopic quantities [45, 48, 46].

Our main finding is the existence of two distinct temperature regions in the AFM phase where Co spins exhibit contrasting low-energy dynamics. In Region III, the spin excitation spectrum has a significant low-energy weight that is enhanced strongly with magnetic field. In contrast, Region IV is likely described by spin-wave excitations. As mentioned above the transition between the two regions is gradual and is possibly a crossover rather than a phase transition with critical dynamics.

A transition to the high-field ordered phase (Regions III^′ and IV^′) was detected around 6 T by the static quantity $B_{\mathrm{hf,pk}}$, reflecting the magnetization jump associated with a reversal of canting moments [45]. The low-energy spin dynamics is essentially unchanged by this transition, although the field response of $1/T_{1}$ is somewhat weakened. A spin-disordered phase does not appear up to 9 T.

The exchange interactions between Co spins have been evaluated by several groups via powder INS techniques, but the results are not settled yet [41, 46, 42, 40]. Here we present an independent evaluation of the dominant interaction strength from the high-$T$ limiting value of ²³Na $1/T_{1}$ which is helpful in justifying a proper energy scale of this compound.

In the PM phase at temperatures much higher than the exchange interactions, the spin dynamics is modeled by Gaussian random modulation of individual spins under the influence of exchange-coupled neighbors. The characteristic (exchange) frequency $\omega_{e}$ is determined solely by the interaction strengths, leading to $T$-independent $1/T_{1}$ [79, 80]. In the presence of both dipolar and isotropic transferred hyperfine coupling, $1/T_{1}$ in the high-$T$ limit is given as a sum of two contributions;

The dipolar contribution is expressed as

where $r_{l}$ is a distance between the nucleus and the $l$-th electron spin, $g$ is the $g$-factor. $\omega_{e}$ is given in terms of the exchange energy $J_{ij}$ between $i$-th and $j$-th electron spins and the number of $j$-th spins $z_{j}$ coupled to the $i$-th spin as

For the contribution of the transferred hyperfine coupling which we assume to come from the nearest-neighbor electron spins, we have

Here $a_{\mathrm{hf}}$ is the hyperfine coupling constant, $z_{\mathrm{n}}$ is the number of nearest-neighbor spins coupled to the nucleus. We take the structural model of Xiao et al. [39] for simplicity and put $z_{\mathrm{n}}=4$.

Evaluating the sum $\sum_{l}r_{l}^{-6}$ in Eq. (3) within a sphere of radius $100\mathrm{~{}\AA}$ and using the value of $a_{\mathrm{hf}}$ determined from the $B_{\mathrm{hf,pk}}$-$M$ scaling in Fig. 4(d), we obtain the ratio of the two contributions as $T_{1\infty,\mathrm{dip}}/T_{1\infty,\mathrm{tr}}=a_{\mathrm{hf}}^{2}/2z_{\mathrm{n}}\sum_{l}r_{l}^{-6}\approx 0.02$. Hence we neglect the contribution $1/T_{1\infty,\mathrm{tr}}$ and go on to evaluate $\omega_{e}$ using Eq. (3). Adopting the isotropic value $g=4.33$ for the $S=j_{\mathrm{eff}}=1/2$ manifold of Co²⁺ for simplicity [20], we finally obtain $\omega_{e}=4.1\times 10^{12}\mathrm{~{}s^{-1}}$ from the observed value $1/T_{1\infty}=100\mathrm{~{}s^{-1}}$.

Since the contributions of $J_{ij}$’s to $\omega_{e}$ are all additive, the maximum value of some $J_{ij}$ may be estimated by neglecting all the other contributions. By putting $J_{ij}=0$ other than the nearest-neighbor Heisenberg coupling $J$ and by putting $z_{1}=3$, we get $\lvert J\rvert=26$ K (2.2 meV). Note that the sign of $J_{ij}$ cannot be determined by the present analysis. In estimating the Kitaev coupling $K$, putting $z_{1}=1$ gives $\lvert K\rvert=44$ K (3.8 meV). The obtained values of $\lvert J\rvert$ and $\lvert K\rvert$ fall in the same order as those evaluated from the INS experiments [41, 42, 40, 46, 47]. According to the microscopic model of the exchange interactions in Na₂Co₂TeO₆ [35], this value of $K$ sets the energy scale $t^{2}/U\approx 1.1$ meV where $t$ is the hopping amplitude and $U$ is Coulomb repulsion defined in Refs. 34, 35.

At temperatures $T_{N}<T\lesssim 2T_{N}$ (Region II), the neutron diffuse scattering experiment has revealed that short-range spin correlations develop within the honeycomb planes [38]. The field suppression of ²³Na $1/T_{1}$ thus suggests a reduction of the in-plane spin correlation length in magnetic fields. In order to quantify the $B$ dependence of $1/T_{1}$ from such a standpoint, we analyze the $T$ dependence of $1/T_{1}$ based on a description of 2D Heisenberg antiferromagnets (HAFs) in terms of the quantum non-linear $\sigma$ model. Similar analysis has recently been applied to ²³Na $1/T_{1}$ in a single crystal of Na₂Co₂TeO₆ at a moderate field [49].

We assume that the system is in the renormalized-classical (RC) regime of 2D HAFs. For the system with collinear order [81, 82, 83], the correlation length $\xi$ grows exponentially with decreasing $T$ in the RC regime as $\xi\propto\exp(2\pi\rho_{s}/T)$ where $\rho_{s}$ is the spin stiffness constant ⁶⁶6Including the correction term of $O(T/2\pi\rho_{s})$ for $\xi$ introduces a minor change of parameters ($\sim$10% reduction of $\rho_{s}$ for example) but does not change the results qualitatively. The spin-lattice relaxation rate is given as $1/T_{1}\propto T^{3/2}\xi$, so that

Considering the possibility of triple-$\mathbf{q}$ order [44], we also examine the non-linear $\sigma$ model developed for the system with noncollinear order [85, 86, 87]. For this type of order, $\xi\propto T^{-1/2}\exp(4\pi\rho_{s}^{\,\prime}/T)$ and $1/T_{1}\propto T^{7/2}\xi$, giving [88]

Here $\rho_{s}^{\,\prime}$ is the corresponding spin stiffness constant.

Figure 8 shows $1/T_{1}T^{3/2}$ plotted against the inverse temperature $1/T$ in a semi-logarithmic scale. $1/T_{1}$ is found to obey the scaling relation Eq. (6) below 50 to 32 K (40 to 24 K at high fields) which ensures that the system is indeed in the RC regime. The scaling form Eq. (7) for noncollinear order is also satisfied in the same $T$ range (not shown). $\rho_{s}$ decreases monotonically with increasing field as shown in the inset of Fig. 8, showing a tendency to vanish at 10$-$12 T. The $B$ dependence of $\rho_{s}^{\,\prime}$ is identical to that of $\rho_{s}$ and is scaled with a multiplicative factor $\rho_{s}/\rho_{s}^{\,\prime}=1.17$ [Fig. 7].

The spin stiffness constant characterizes the rigidity of an ordered state and is nonzero only in a phase with magnetic LRO. In the quantum non-linear $\sigma$ model for 2D HAFs, it is renormalized by quantum fluctuations and vanishes on approaching a quantum critical point, beyond which a quantum disordered state appears [81]. The field-induced reduction of $\rho_{s}$ and $\rho_{s}^{\,\prime}$ thus suggests the system getting closer to the quantum disordered phase present above a certain critical field. This also implies that the high-field phase of Na₂Co₂TeO₆ is not a QSL but a partially polarized phase showing bosonic excitations. As a rough estimate of the critical field $B_{c}$, we fitted the $B$ dependence of $\rho_{s}$ and $\rho_{s}^{\,\prime}$ to the power law $\rho_{s},\rho_{s}^{\,\prime}\propto(B_{c}-B)^{p}$, getting $B_{c}=10.4$ T and $p=0.22$. The obtained value of $B_{c}$ agrees well with the value $B_{c}\approx 10$ T determined from the magnetization measurement [48], although the agreement seems fortuitous considering the lack of our data closer to $B_{c}$.

The $B$ dependence of $\rho_{s}$ and $\rho_{s}^{\,\prime}$ is also shown in Fig. 7 for comparison with the phase diagram. Surprisingly enough, $\rho_{s}$ and $\rho_{s}^{\,\prime}$ trace $T_{N}$ as a function of magnetic field with appropriate scale factors. ($\rho_{s}\approx 0.5T_{N}$.) This implies that the 2D spin correlations and the 3D magnetic LRO are characterized by a common energy scale that is renormalized by magnetic field. It is plausible that the magnetic field changes a balance of competing interactions to make the system more frustrated and suppress magnetic LRO, which may appear as a decrease of the characteristic energy and the spin stiffness constant, leading to the suppression of $\xi$ and $1/T_{1}$. In fact, frustration reduces the spin stiffness constant in 2D quantum HAFs [89, 90]. It is, however, puzzling that the macroscopic susceptibility scarcely depends on the magnetic field in the PM phase.

$1/T_{1}$ deviate from the RC scaling relations Eqs. (6) and (7) close to $T_{N}$, signaling an onset of 3D critical slowing down. One may expect a power-law dependence of $1/T_{1}$ on the reduced temperature $\epsilon=T/T_{N}-1$ in the critical region, the exponent of which provides information on the universality class of the phase transition. We do not pursuit this subject because ²³Na $1/T_{1}$ is strongly modified by applying field in the critical region, which prevents us from extracting a reliable value of the critical exponent.

One of the most important characteristics of the spin dynamics in Na₂Co₂TeO₆ is the existence of an intermediate temperature region in the AFM phase (Region III) in which Co spins exhibit unconventional dynamics. The magnetic excitation spectrum in Region III ($T_{N}/2\lesssim T<T_{N}$) comprises of a broad continuum centered at the ordering wave vector $\mathbf{Q}=(\frac{1}{2},0,0)$ rather than a distinct spin-wave mode [44]. At low energies probed via $1/T_{1}$, there remains an unusually large spectral weight of spin fluctuations enhanced strongly with magnetic field. In contrast, the spin dynamics at low $T$ in Region IV ($T\lesssim T_{N}/2$) looks more conventional. The continuum is replaced with a gapped spin-wave mode, and $1/T_{1}$ exhibits a rapid decrease on cooling which evidences the disappearance of low-energy spectral weight. The key ingredients in understanding the spin dynamics of Na₂Co₂TeO₆ would thus be (i) the origin of a continuum and the large spectral weight at low energies in Region III, (ii) the origin of field enhancement of the low-energy spectral weight, and (iii) the trigger for the formation of a spin-wave mode in Region IV.

The most likely cause of the broad continuum accompanying a large spectral weight at low energies would be strong damping of spin waves. The appearance of long-lived spin waves at low $T$ is then understood as resulting from diminished damping. It is apparent that conventional mechanisms for spin-wave damping in collinear antiferromagnets do not apply because they yield too small damping to account for a broad feature of the excitation spectrum in Region III [91]. If allowed by symmetry, cubic anharmonicities in the magnon effective Hamiltonian enable the coupling between one- and two-magnon states which possibly leads to strong damping of the one-magnon mode degenerate with a two-magnon continuum [92, 93]. The cubic terms exist in noncollinear antiferromagnets as well as the Kitaev magnets with off-diagonal interactions and may lead to severe damping at relatively high energies, leaving a low-energy one-magnon mode almost untouched. This is not the case in Na₂Co₂TeO₆ because there is no distinct excitation mode at low energies in Region III on the one hand, and the damping is not so severe at high energies in Region IV on the other [44]. Among other things, the fact that the damping changes rapidly its character around a certain temperature ($\sim T_{N}/2$) seems difficult to be accounted for by a known mechanism for damping which usually gives a smooth variation of magnon lifetime with $T$.

Although a prime mechanism for the spin-wave damping is unidentified, temperature and field variation of the damping may be argued qualitatively based on the general expression of nuclear spin-lattice relaxation rate [79],

Here $\chi^{\prime\prime}(\mathbf{q},\omega)$ is the imaginary part of the dynamical spin susceptibility, $\omega_{0}$ is the nuclear Larmor frequency, and $A(\mathbf{q})$ is the hyperfine form factor determined by a geometry of the nuclear site. Equation (8) tells us that $1/T_{1}$ is determined by a spectral weight of spin fluctuations at a very low frequency $\omega_{0}\sim 10^{8}\mathrm{~{}s^{-1}}$ ($\hbar\omega_{0}\sim 0.1\mathrm{~{}\mu eV}$). Since there appears a staggered hyperfine field at the Na sites below $T_{N}$ as revealed from the NMR line broadening, $A(\mathbf{Q})$ at the ordering wave vector $\mathbf{Q}$ is nonzero and a dominant contribution to $1/T_{1}T$ will come from $\mathbf{q}\sim\mathbf{Q}$ in the AFM phase.

The most intuitive view of the results not relying on the specific model is to interpret $1/T_{1}T$ as a measure of low-energy spin excitations represented by $\chi^{\prime\prime}(\mathbf{q},\omega_{0})$. It is thus apparent from Figs. 6(b) and 7 that the low-energy spin excitations remain enhanced in Region III in the AFM phase, especially at high fields of 6$-$7 T as pointed out in the preceding section. It is also obvious that the active excitation channels survive to lower temperatures at higher fields. The excitation channels activated in Region III almost disappear on entering Region IV as evidenced by a steep decrease of $1/T_{1}T$ below $\sim$10 K.

If we take a model of damped harmonic oscillator for spin-wave excitations [94], we may obtain semi-quantitative information on the spin-wave damping. According to the model, we have a contribution of the dynamical susceptibility to $1/T_{1}T$ in the limit of $\omega_{0}\rightarrow 0$ as $\chi^{\prime\prime}(\mathbf{q},\omega_{0})/\omega_{0}\propto\gamma_{\mathbf{q}}/\omega_{\mathbf{q}}^{2}$ where $\omega_{\mathbf{q}}$ is an undamped spin-wave frequency and $\gamma_{\mathbf{q}}$ is a damping constant. $\omega_{\mathbf{q}}$’s are expected to depend only weakly on $T$ except in the vicinity of $T_{N}$ and unless there is a drastic change of the magnetic structure. The $T$ dependence of $1/T_{1}T$ at a constant $B$ is thus dominated by that of $\gamma_{\mathbf{q}}$. As to the $B$ dependence, both $\omega_{\mathbf{q}}$ and $\gamma_{\mathbf{q}}$ may vary with $B$ and affect the behavior of $1/T_{1}T$ because of unknown effects of magnetic field on them.

As shown in Figs. 6(b) and 7, $1/T_{1}T$ is peaked around 6$-$7 T at a constant $T$ in Region III. This suggests an increase of $\gamma_{\mathbf{q}}$ and/or a decrease of $\omega_{\mathbf{q}}$ with magnetic field. Notice that the latter matches the field-induced reduction of the characteristic energy inferred from the RC analysis of $1/T_{1}$ in the PM phase. At a constant $B$, on the other hand, the $T$ dependence of $1/T_{1}T$ should be ascribed to that of $\gamma_{\mathbf{q}}$ as mentioned above. Therefore, a broad hump of $1/T_{1}T$ suggests enhancement of the spin-wave damping at high fields toward the boundary between Regions III and IV. The origin of such strong damping and its unusual $T$ and $B$ dependence is unclear. Frustration might play an important role in making the spin excitations incoherent as observed in the triangular-lattice antiferromagnet NaCrO₂ below the spin-freezing temperature $T_{c}\approx 41$ K [95]. Note, however, that the excitation spectrum of NaCrO₂ becomes dispersive below about $0.75T_{c}$ triggered by the onset of short-range 3D spin correlations, whereas in Na₂Co₂TeO₆ both in-plane and out-of-plane spin correlations show no appreciable change across the boundary between Regions III and IV [38].

The presence of magnetic scattering centers is another possibility of strong spin-wave damping. This may provide a reasonable account for the broad continuum in Region III that changes to the dispersive mode in Region IV as described below. Structural disorder in Na layers and a related distribution of the interlayer magnetic coupling seem less important because they would give $T$-independent damping. It is worth mentioning here that the qualitative change of the excitation spectrum with $T$ in Na₂Co₂TeO₆ resembles the behavior observed in some geometrically frustrated magnets such as the kagomé staircase Ni₃V₂O₈ [96] and the pyrochlore Gd₂Ti₂O₇ [97]. These compounds have an intermediate temperature phase just below $T_{N}$ in which the spin excitations are quasielastic or show only broad features and a low temperature phase with collective excitations. The intermediate phases are identified as a partially-disordered state where long-range ordered moments coexist with disordered (paramagnetic) moments, whereas all the sites are ordered in the low temperature phase. This suggests a vital role of partially-disordered moments in the spin-wave damping in these materials. Such a scenario may be realized in Na₂Co₂TeO₆ if the AFM phase has triple-$\mathbf{q}$ rather than single-$\mathbf{q}$ zigzag order [44].

In the triple-$\mathbf{q}$ ordered state shown in Fig. 1(b), three quarters of the Co atoms show noncollinear spin arrangement with a vortex-like texture in the honeycomb planes, while the remaining Co atoms become “spinless”, which means that the ordered moment is absent, or have only an out-of-plane Néel component. The mean fields originating from the in-plane component of the majority spins cancel at the minority site, which may allow the minority spin to fluctuate in a large amplitude. The absence of the signal wipeout implies that if the minority spins are not ordered in Region III, they fluctuate in a time scale much faster than the NMR time scale like a paramagnetic moment. Spin waves propagating on the majority sites would be strongly damped to give a broad continuum if the resulting quasielastic mode of the minority spins overlaps energetically with the spin-wave mode.

Since the minority spins are coupled via the third-neighbor interactions and possibly via the effective interactions mediated by the majority spins due to quantum fluctuations around the order [98, 99], they would participate in magnetic LRO at low enough temperatures exhibiting slowing of spin fluctuations. A broad hump of the magnetic specific heat observed around 10 K [45, 46] might be associated with such a change of the minority spin state. The spin-wave damping would be diminished as a spectral weight of the quasielastic mode shifts to lower energies, restoring collective excitations at low $T$.

The low-energy spectral weight of the majority spins arises from the spin-wave damping and is reduced to give a decreasing contribution to $1/T_{1}T$ at low $T$. On the other hand, slowing of the minority spin fluctuations would contribute a peak or hump of $1/T_{1}T$ like the case of a magnetic phase transition. The $T$ and $B$ dependence of $1/T_{1}T$ should be determined by a balance between the two contributions. Actually, a feeble anomaly of $1/T_{1}T$ around 15 K at low fields implies the crossover nature of slowing down rather than a sharp transition with critical dynamics as noted in the previous section. A gradual increase of the NMR line splitting below $\sim$15 K [49] might be related to this crossover and the resulting appearance of a static moment on the minority site. The broad hump of $1/T_{1}T$ around 10 K at high fields should then be ascribed primary to the majority spins suffering strong damping even at that temperature for unknown reasons, although it might be possible that the minority spin dynamics becomes more critical at higher fields to contribute the hump. The $T$ and $B$ dependence of $1/T_{1}T$ around the boundary between Regions III and IV is complex and is not fully understood in terms of the minority spin ordering. Further investigations are needed to clarify the field-dependent spin dynamics in the AFM phase of Na₂Co₂TeO₆.

It has recently been argued whether Na₂Co₂TeO₆ serves as a canonical example of the Kitaev magnet. As described in the preceding sections, the spin dynamics in Na₂Co₂TeO₆ displays distinct features from those in other Kitaev candidates such as $\alpha$-RuCl₃ and Na₂IrO₃. The magnetic excitation spectrum in the AFM phase of Na₂Co₂TeO₆ is characterized by the sole existence of a broad continuum or a distinct spin-wave mode, both of which has a dominant intensity around the M points of the 2D Brillouin zone ($\mathbf{Q}=(\frac{1}{2},0,0)$ and the equivalents) [44]. On the other hand, a continuum in $\alpha$-RuCl₃ and Na₂IrO₃ is centered at the $\Gamma$ point (zone center) and coexists with spin-wave modes below $T_{N}$ [23, 25]. The major excitations around the M point in Na₂Co₂TeO₆ are possibly ascribed to large third-neighbor coupling $J_{3}$ suggested from the powder INS [40, 46, 42, 47]. This is known to stabilize zigzag order but to counteract the formation of a Kitaev QSL [34, 35]. In the context of Kitaev physics, the honeycomb cobaltate Na₃Co₂SbO₆ seems more promising because it exhibits intense excitations around the $\Gamma$ point probably due to smaller $J_{3}$ [41, 42].

On the field evolution of the low-energy spin dynamics, it is interesting to compare our results of $1/T_{1}$ with those for other Kitaev candidates exhibiting a similar field-induced transition from an AFM phase to a spin-disordered phase. To the best of the authors’ knowledge, $\alpha$-RuCl₃ is only one such example with extensive field-dependent NMR studies [55, 56, 57, 59]. The AFM phase of another well-studied candidate Na₂IrO₃ is robust against a field [43]. In fact, $1/T_{1}$ at the Na sites in Na₂IrO₃ is insensitive to field and shows conventional behaviors above and below $T_{N}$ [58].

$\alpha$-RuCl₃ seems to be the best reference as it shows zigzag order like Na₂Co₂TeO₆. The in-plane critical field $B_{c}$ to the high-field disordered phase is 7$-$8 T [28, 29, 30]. Most of the NMR experiments on $\alpha$-RuCl₃, however, focused on the behavior near and above $B_{c}$ and the field evolution of $1/T_{1}$ has not been investigated systematically in the low-field region. Although a direct comparison of the results is limited to a narrow field range, the field response of $1/T_{1}$ has distinct differences between the two compounds. This may reflect the presence of the intermediate temperature region (Region III) characterized by a substantial low-energy spectral weight in Na₂Co₂TeO₆ and the excitation continuum coexisting with spin-wave modes in $\alpha$-RuCl₃. Indeed, $1/T_{1}$ at the ³⁵Cl site in the AFM phase of $\alpha$-RuCl₃ is relatively insensitive to field except above the field $B_{c}^{\prime}=7.1$ T $(<B_{c})$ where gapless magnon excitations have been suggested [59]. This shows a marked contrast to strong field enhancement of ²³Na $1/T_{1}$ in Region III of Na₂Co₂TeO₆ starting far below the critical field $B_{c}\approx 10$ T. The field-insensitive response of $1/T_{1}$ at $B<B_{c}^{\prime}$ in $\alpha$-RuCl₃ would be due to gapped magnon excitations, but at $T<4$ K well below $T_{N}=6.5$ K, $1/T_{1}$ is contributed by a residual mode that grows with field on crossing $B_{c}^{\prime}$ and becomes dominant above $B_{c}$ [55, 59]. Such a mode was not detected in Na₂Co₂TeO₆ down to 3.5 K and may be associated with the continuum around the $\Gamma$ point in $\alpha$-RuCl₃ identified as excitations inherent to Kitaev QSLs.

Despite the apparent differences in the low-energy spin dynamics probed by $1/T_{1}$, the two compounds have a lot of similarities in their response to magnetic field; closing of the AFM phase suggesting the existence of a quantum critical point, the possible appearance of a field-induced QSL phase, and so on. It is often encountered in strongly frustrated magnets that the magnetic LRO is controlled by sub-leading interactions instead of the leading one like the Kitaev coupling and by external perturbations such as magnetic field, pressure, and in some cases spin defects. Since the low-energy sector of a magnetic excitation spectrum is very much affected and reconstructed by these interactions and perturbations, there will be a wide variety of phases and behaviors in real materials which at first glance look very different. Concerning the present case, it might be possible that the unusual temperature and field evolution of low-energy spin dynamics in Na₂Co₂TeO₆ is described in terms of a generalized Kitaev model by including necessary factors. It is worth noting that for the reported values of the Kitaev coupling $K$ [41, 42, 40, 47], Region III is around or lower than the crossover temperature $T_{H}\sim 0.375K$ below which localized $Z_{2}$ fluxes and itinerant Majorana fermions are expected to emerge [100]. The close resemblance between $1/T_{1}T$ and $C_{m}/T$, however, implies confinement of the fractionalized particles to magnons. Our findings on the low-energy spin dynamics of Na₂Co₂TeO₆ will thus provide new insights into Kitaev-derived spin models as well as more conventional models on the honeycomb lattice, promoting future studies in this research field. From an experimental side, field-dependent microscopic measurements complementing NMR, such as neutron, Raman, and terahertz spectroscopies, are highly required.