Majorana dynamical mean-field study of spin dynamics at finite temperatures in the honeycomb Kitaev model

A prominent feature of quantum spin liquids is fractionalization of the spin degree of freedom. Fractionalized excitations have their own dynamics in different energy scales, and hence, affect finite-temperature ($T$) properties in a peculiar manner even in the paramagnetic state harboring the quantum spin liquid state. We here present a comprehensive theoretical study of the spin dynamics in a wide $T$ range for the Kitaev model on a honeycomb lattice, whose ground state is a quantum spin liquid. In this model, the fractionalization occurs to break up quantum spins into itinerant matter fermions and localized $Z_2$ fluxes, which results in two crossovers at very different $T$ scales. Extending the previous study for the isotropic coupling case [J. Yoshitake, J. Nasu, and Y. Motome, Phys. Rev. Lett. 117, 157203 (2016)], we calculate the dynamical spin structure factor $S(\mathbf{q},\omega )$, the NMR relaxation rate $1/T_1$, and the magnetic susceptibility $\chi$ while changing the anisotropy in the exchange coupling constants, by using the dynamical mean-field theory based on a Majorana fermion representation. We describe the details of the methodology including the continuous-time quantum Monte Carlo method for computing dynamical spin correlations and the maximum entropy method for analytic continuation. We confirm that the combined method provides accurate results in a wide $T$ range including the region where the spins are fractionalized. We find that also in the anisotropic cases the system exhibits peculiar behaviors below the high-$T$ crossover whose temperature is comparable to the average of the exchange constants: $S(\mathbf{q},\omega )$ shows an inelastic response at the energy scale of the averaged exchange constant, $1/T_1$ continues to grow even though the equal-time spin correlations are saturated and almost $T$ independent, and $\chi$ deviates from the Curie-Weiss behavior. In particular, when the exchange interaction in one direction is stronger than the other two, the dynamical quantities exhibit qualitatively different $T$ dependences from the isotropic case at low $T$, reflecting the opposite parity between the flux-free ground state and the flux-excited state, and a larger energy cost for flipping a spin in the strong interaction direction. On the other hand, when the exchange anisotropy is in the opposite way, the results are qualitatively similar to those in the isotropic case. All these behaviors manifest the spin fractionalization in the paramagnetic region. Among them, the dichotomy between the static and dynamical spin correlations is unusual behavior hardly seen in conventional magnets. We discuss the relation between the dichotomy and the spatial configuration of the $Z_2$ gauge fields. Our results could stimulate further experimental and theoretical analyses of candidate materials for the Kitaev quantum spin liquids.

Figure 1

(a) Schematic picture of the Kitaev model on the honeycomb lattice. The blue, green, and red bonds represent the $p=x,y$, and $z$ bonds in Eq. (1), respectively. The dashed oval represents the 26-site cluster used in the CDMFT calculations. (b) The first Brillouin zone (black hexagon) and the symmetric lines (red lines) used in Figs. 3 and 4.

Figure 2

The specific heat $C_v$ and equal-time spin correlations for the NN sites, ${\langle S_j^z{S}_{j^{\prime }}^z\rangle }_{\mathrm{NN}}$ and ${\langle S_j^x{S}_{j^{\prime }}^x\rangle }_{\mathrm{NN}}$, obtained by the Majorana CDMFT for the FM case at (a) $\alpha =0.8$ and (b) 1.2. Note that ${\langle S_j^x{S}_{j^{\prime }}^x\rangle }_{\mathrm{NN}}={\langle S_j^y{S}_{j^{\prime }}^y\rangle }_{\mathrm{NN}}$ from the symmetry. QMC data in Ref. [10] are plotted by gray symbols for comparison.

Figure 3

Dynamical spin structure factor $S(\mathbf{q},\omega )$ obtained by the Majorana CDMFT+CTQMC method for the FM case with (a), (d), (g), and (j) $\alpha =1.0$; (b), (e), (h), and (k) $\alpha =0.8$; and (c), (f), (i), and (l) $\alpha =1.2$; as well as (a), (b), and (c) $T\simeq 2T_{\mathrm{L}}$; (d), (e), and (f) $T\simeq \sqrt{T_{\mathrm{L}}{T}_{\mathrm{H}}}$; (g), (h), and (i) $T\simeq 0.64T_{\mathrm{H}}$; and (j), (k), and (l) $T\simeq 6.4T_{\mathrm{H}}$. Here, $T_{\mathrm{L}}\simeq 0.012$, 0.0052, and 0.0075 for $\alpha =1.0$, 0.8, and 1.2, respectively, while $T_{\mathrm{H}}\simeq 0.375$ for all the cases.

Figure 4

Dynamical spin structure factor $S(\mathbf{q},\omega )$ obtained by the Majorana CDMFT+CTQMC method for the AFM case. The values of $\alpha$ and $T$ are common to Fig. 3.

Figure 5

(a) $S(\mathrm{\Gamma },\omega )$, (c) $S({\mathrm{K}}_1,\omega )$, and (e) $S({\mathrm{K}}_2,\omega )$ for the FM case with $\alpha =0.8$ at several $T$. The corresponding contour plots in the $T\text{−}\omega$ plane are shown in (b), (d), and (f). The arrows indicate the temperatures used for the data in (a), (c), and (e), while the white and gray dotted lines indicate $T_{\mathrm{H}}$ and $T_{\mathrm{L}}$, respectively. Note that the $T$ set is common to that used in Figs. 3 and 4. The dashed curve in (b) represent the average frequency of $S(\mathrm{\Gamma },\omega )$ (see the text for details). In (a), (c), and (e), the error bars are shown for every ten data along the $\omega$ axis.

Figure 6

(a) $S(\mathrm{\Gamma },\omega )$, (c) $S({\mathrm{K}}_1,\omega )$, and (e) $S({\mathrm{K}}_2,\omega )$ for the AFM case with $\alpha =0.8$ at several $T$. The corresponding contour plots in the $T\text{−}\omega$ plane are shown in (b), (d), and (f). The notations are common to those in Fig. 5.

Figure 7

(a) $S(\mathrm{\Gamma },\omega )$, (c) $S({\mathrm{K}}_1,\omega )$, and (e) $S({\mathrm{K}}_2,\omega )$ for the FM case with $\alpha =1.2$ at several $T$. The corresponding contour plots in the $T\text{−}\omega$ plane are shown in (b), (d), and (f). The notations are common to those in Fig. 5.

Figure 8

(a) $S(\mathrm{\Gamma },\omega )$, (c) $S({\mathrm{K}}_1,\omega )$, and (e) $S({\mathrm{K}}_2,\omega )$ for the AFM case with $\alpha =1.2$ at several $T$. The corresponding contour plots in the $T\text{−}\omega$ plane are shown in (b), (d), and (f). The notations are common to those in Fig. 5.

Figure 9

$T$ dependences of the NMR relaxation rate $1/T_1^p$ ($p=z,x$) at (a) $\alpha =1.0$, (b) 0.8, and (c) 1.2. Note that $1/T_1^z=1/T_1^x$ for $\alpha =1.0$ and $1/T_1^x=1/T_1^y$ for all the cases from the symmetry. The vertical dotted lines indicate $T_{\mathrm{L}}$ and $T_{\mathrm{H}}$ for each $\alpha$.

Figure 10

$T$ dependences of the magnetic susceptibility ${\chi }^p$ ($p=z,x$) at (a) $\alpha =1.0$, (b) 0.8, and (c) 1.2 for the FM case. The dashed curves represent ${\chi }_{\mathrm{CW}}^p$ in Eq. (45). The red and black dashed-dotted curves in (b) represent ${\chi }_{\mathrm{dimer}}^p$ for $p=z$ [Eq. (46)] and $p=x$ [Eq. (47)], respectively. Note that ${\chi }^z={\chi }^x$ for $\alpha =1.0$ and ${\chi }^x={\chi }^y$ for all the cases from the symmetry. The vertical dotted lines indicate $T_{\mathrm{L}}$ and $T_{\mathrm{H}}$ for each $\alpha$.

Figure 11

$T$ dependences of the magnetic susceptibility ${\chi }^p$ ($p=z,x$) at (a) $\alpha =1.0$, (b) 0.8, and (c) 1.2 for the AFM case. The notations are common to those in Fig. 10.

Figure 12

(a)–(c) $4{\langle S_j^p{S}_{j^{\prime }}^p\rangle }_{\mathrm{NN}}$ for the FM case, (d)–(f) the on-site components of $1/T_1^p$, and ${\chi }^p$ for the (g)–(i) FM and (j)–(l) AFM cases ($p=z,x$) calculated by setting all ${\eta }_r=1$ (uniform) and all ${\eta }_r$ being random (random) in the CDMFT calculations: (a), (d), (g), and (j) $\alpha =1.0$; (b), (e), (h), and (k) $\alpha =0.8$; and (c), (f), (i), and (l) $\alpha =1.2$. The data in (h) are multiplied by 0.25 and 3 for $p=z$ and $x$, respectively. For $\alpha =1.0$, the results are equivalent for $p=x,z$. For comparison, we plot the data in Figs. 2 (CDMFT), 9, 10, and 11 (CDMFT+CTQMC). The vertical dotted lines represent $T_{\mathrm{L}}$ and $T_{\mathrm{H}}$ for each $\alpha$.

Figure 13

Correlation functions between the $Z_2$ gauge fields on different $z$ bonds: (a) $\alpha =1.0$, (b) 0.8, and (c) 1.2. The correlations are measured for ${\eta }_r$ and ${\eta }_{r^{\prime }}$ which are sandwiched by the common two $xy$ chains; $|r-r^{\prime }|$ is the distance in the unit of the width of the hexagon [see the inset in (a) for the example of $|r-r^{\prime }|=3]$. The yellow hexagons represent two neighboring fluxes which are flipped by an operation of $S_j^z$ on the sharing $z$ bond; see the text for the details. The vertical dotted lines represent $T_{\mathrm{L}}$ and $T_{\mathrm{H}}$ for each $\alpha$.

Figure 14

Schematic pictures of the different types of clusters used in the benchmark of CDMFT. The color of the bonds are common to Fig. 1.

Figure 15

Cluster-size dependences of the magnetic susceptibility ${\chi }^p$ for (a)–(e) the FM case and (f)–(j) the AFM case, and (k)–(o) the on-site component of the NMR relaxation rate $1/T_1^p$: (a), (f), and (k) $\alpha =1.0$; (b), (c), (g), (h), (l), and (m) $\alpha =0.8$; and (d), (e), (i), (j), (n), and (o) $\alpha =1.2$. The data for two different $T$ are plotted in each case. In (a), (f), and (k), the data are common to $p=z$ and $x$. Calculations are performed for the cluster series denoted in (A) Fig. 1, (B) Fig. 14, and (C) Fig. 14. Symbols in (a)–(e) are common for the same parameters in (f)–(o).

Figure 16

Comparison between the MEM, RTE, and exact-DOS results for (a) and (b) $S_{j,j}^x\left(\omega \right)$; (c) and (d) $S_{\mathrm{NN}}^x\left(\omega \right)$; (e) and (f) $S_{j,j}^x\left(\omega \right)-S_{\mathrm{NN}}^x\left(\omega \right)$; and (g) and (h) $S_{j,j}^z\left(\omega \right)$ at (a), (c), (e), and (g) $T=0.375$ and (b), (d), (f), and (h) $T=0.0375$.

Figure 17

$S^z(\mathrm{\Gamma },\omega )$ obtained by the Majorana CDMFT+CTQMC method for the FM case for $\alpha =0.8$ at $T=0.003$. Vertical line at $\omega \sim 0.042$ represents the value of flux gap in the ground state.

Figure 18

Spin correlations as functions of $T$ and $\omega$ for the FM case: (a), (c), (e), (g), and (i) on-site components $S_{j,j}\left(\omega \right)$ and (b), (d), (f), (h), and (j) NN-site components $S_{\mathrm{NN}}\left(\omega \right)$. (a) and (b) are for $S_{j,j}^x\left(\omega \right)=S_{j,j}^z\left(\omega \right)$ and $S_{\mathrm{NN}}^x\left(\omega \right)=S_{\mathrm{NN}}^z\left(\omega \right)$, respectively, at $\alpha =1.0$. (c), (d), (e), and (f) [(g), (h), (i), and (j)] are for $S_{j,j}^x\left(\omega \right)$, $S_{\mathrm{NN}}^z\left(\omega \right)$, $S_{j,j}^z\left(\omega \right)$, and $S_{\mathrm{NN}}^z\left(\omega \right)$, respectively, at $\alpha =0.8$ (1.2). The white and gray dotted lines indicate $T_{\mathrm{H}}$ and $T_{\mathrm{L}}$, respectively, for each $\alpha$.

Figure 19

$T$ dependences of the Korringa ratio ${\mathcal{K}}^p=1/\left(T_1^{p}T{\left({\chi }^p\right)}^2\right)$ for the FM case at (a) $\alpha =1.0$, (b) 0.8, and (c) 1.2 ($p=z,x$). Note that ${\mathcal{K}}^z={\mathcal{K}}^x$ for $\alpha =1.0$ and ${\mathcal{K}}^x={\mathcal{K}}^y$ for all the cases from the symmetry. The vertical dotted lines indicates $T_{\mathrm{L}}$ and $T_{\mathrm{H}}$ for each $\alpha$.

Figure 20

$T$ dependences of the Korringa ratio ${\mathcal{K}}^p=1/\left(T_1^{p}T{\left({\chi }^p\right)}^2\right)$ for the AFM case at (a) $\alpha =1.0$, (b) 0.8, and (c) 1.2 ($p=z,x$). The notations are common to those in Fig. 19.