Phonon dynamics in the Kitaev spin liquid

The search for fractionalization in quantum spin liquids largely relies on their decoupling with the environment. However, the spin-lattice interaction is inevitable in a real setting. While the Majorana fermion evades a strong decay due to the gradient form of spin-lattice coupling, the study of the phonon dynamics may serve as an indirect probe of fractionalization of spin degrees of freedom. Here we propose that the signatures of fractionalization can be seen in the sound attenuation and the Hall viscosity. Despite the fact that both quantities can be related to the imaginary part of the phonon self-energy, their origins are quite different, and the time-reversal symmetry breaking is required for the Hall viscosity. First, we compute the sound attenuation due to a phonon scattering off of a pair of Majorana fermions and show that it is linear in temperature ($\sim T$). We argue that it has a particular angular dependence providing the information about the spin-lattice coupling and the low-energy Majorana fermion spectrum. The observable effects in the absence of time-reversal symmetry are then analyzed. We obtain the phonon Hall viscosity term from the microscopic Hamiltonian with time-reversal symmetry breaking term. Importantly, the Hall viscosity term mixes the longitudinal and transverse phonon modes and renormalize the spectrum in a unique way, which may be probed in spectroscopy measurement.


I. INTRODUCTION
Quantum spin liquids (QSLs), a particularly fascinating class of frustrated magnets, have been a focus of condensed matter research since the initial proposal [1].These systems evade magnetic order down to zero temperature and harbor a remarkable set of collective phenomena, including topological ground-state degeneracy, long-range entanglement, and fractionalized excitations [2][3][4][5].Of particular interest is the Kitaev honeycomb model which describes a system of spin-1/2 at sites of a honeycomb lattice interacting via Ising-like frustrated nearest-neighbor exchange interactions [3].This model is not only exactly solvable with a QSL ground state, but is also realizable in real materials [6][7][8][9][10][11][12].
Here we propose to study the signatures of the fractionalized excitations in Kitaev materials by exploiting their coupling to lattice vibrations.Indeed, it is well known that the coupling between collective spin and lattice degrees of freedom plays an important role on the fundamental properties of correlated materials, and in many cases, it has been explored rather carefully.For example, attenuation of sound due to electron-phonon coupling has been widely used to study the electronic properties and phase transitions of various solids with complex order parameters [30][31][32] and was proved to provide useful information on both normal and super-conducting states [33][34][35].The magneto-elastic coupling also plays an important role in probing vestigial phases of frustrated magnets [36].The possibility to use the sound attenuation to study the nature of two-dimensional frustrated magnetic systems, and in particular QSLs, have been also recently discussed in the literature [37][38][39][40].In particular, the spin-lattice coupling was shown to play an essential role in the theoretical interpretation of some experiments in the study of Kitaev materials, e.g. the thermal Hall transport measurement [41][42][43].This suggests that, in a realistic setting, the spin-lattice coupling in Kitaev materials could be large enough to affect the phonon dynamics [40].However, a detailed analysis of the sound attenuation in Kitaev magnets, and in particular the signatures constrained by the symmetry of the spin-lattice coupled system, remain mostly unexplored.
In this paper we argue that the phonon dynamics can be used to probe spin fractionalization in Kitaev materials, and in particular in α-RuCl 3 [18,41,[44][45][46].As a proof of principle, we study the pure Kitaev model [3] in the isotropic limit on a 2D honeycomb lattice [see Fig. 1 (a)], with parameters such as sound velocity v s and Kitaev interaction J K extracted from experimental data on α-RuCl 3 [8,46].In the pure Kitaev model, the spins fractionalize into two types of elementary excitations -Majorana fermions and emergent gauge fluxes [3].The ground state of the isotropic Kitaev spin liquid corresponds to a fixed zero-flux configuration.Therefore at temperatures below the flux gap ∆ flux , the low-energy magnetic excitations are solely dispersive Majorana fermions.We assume that the Majorana fermion-phonon coupling arises from the fact that the Kitaev interaction J K depends on the relative positions between the spins.In order to obtain an effective low-energy theory, we performed a microscopic analysis of the change of the spin exchange energy due to the lattice distortion and obtain the explicit Ω ⊥ H (q) Ω ⊥ 0 (q) Ω  H (q) Ω  0 (q) δ c v ⊥ q/J K 2 ), and Mz = (0, −1) connect nearest neighbors by x, y, and z bonds, respectively.The three sites r, r , r form a representative three-spin link r, r , r yx as defined below Eq. ( 2).(b) Angular dependence of the sound attenuation coefficient α s for the longitudinal mode.α ⊥ s for the transverse mode follows the same shape but is rotated by π/6.(c) Acoustic phonon spectrum Ω(q).The dashed lines depicts the bare spectrum with longitudinal and transverse modes.Solid lines depict the spectrum renormalized due to the mixing of transverse and longitudinal modes when time-reversal symmetry is broken.δc denotes the characteristic dimensionless scale beyond which the phonon spectrum is significantly bended [see Eq. ( 61)].
form of the Majorana fermion-phonon coupling by considering acoustic phonon modes coupled to low-energy spin degrees of freedom expressed in terms of Majorana fermions.We also found that in the low-energy limit, this coupling has essentially the same form as that obtained from the symmetry considerations in Ref. [39] for algebraic spin liquids.
In this picture, we study the observable consequences of the spin-lattice coupling through the phonon dynamics.By calculating the phonon self-energy, we compute the sound attenuation coefficient and the Hall viscosity, analyzing their observable signatures.
We first show that the sound attenuation is determined by the decay of a phonon into a pair of Majorana fermions [see Fig. 2(a)], with the attenuation rate linear in temperature due to the vanishing density of states at the Dirac points.Importantly, this is the dominant process compared with the sound attenuation due to, e.g.phononphonon interactions that scales as ∝ T 4 [47].Moreover, we find that due to the anisotropic form of the Majorana fermion-phonon coupling and the Dirac fermion like low-energy Hamiltonian of the magnetic excitations, the sound attenuation shows a strong angular dependence at the leading order in phonon momentum q [see Fig. 1(b)].Consequently, the angular dependence may offer a quite powerful probe of Majorana fermions with Dirac spectrum in Kitaev materials.
The same Majorana fermion-phonon interaction also gives rise to the finite life-time of the Majorana fermions.This effect, however, is quite weak, and the life-time scales as τ f ∼ Im Σ mf ∼ T 2 , which is much smaller than the typical fermion energy ∼ T .
We next study the modifications of the phonon dynamics in the absence of time-reversal symmetry due to applying a magnetic field h.Assuming the phonons do not couple to the magnetic field directly, these corrections are induced by the spin-lattice coupling.While the sound attenuation coefficient α s is not changed qualitatively in a small magnetic field h 3 /J 2 K v F q/2, the Majorana fermion-phonon coupling introduces another interesting effect to the phonon system, i.e. a Berry phase term that mixes the transverse and longitudinal phonon modes.It is encoded in the phonon long-wavelength effective action as the Hall viscosity term, which is the leading order term breaking the time-reversal symmetry [48,49].It is associated with the non-dissipative response to a velocity gradient of a fluid [48][49][50], and also contributes to the thermal Hall effect [43,51].As it comes from off-shell processes, both the high-energy and low-energy Majorana fermion modes contribute.In our analysis presented below, we calculate the Hall coefficient from the whole Majorana fermion spectrum.We compute the nonperturbative contribution in h to the Hall viscosity coefficient, which matches the result in an infinitesimally small field.We also compute the perturbative correction and show that the Hall coefficient decreases as the magnetic field strength increases.We also note that the mixing renormalizes the phonon spectrum prominently above a characteristic phonon momentum q c [see Fig. 1(c)] that may be observed in some spectroscopy measurements.
The rest of the paper is organized as follows.In Sec.II, we present the spin-phonon Hamiltonian.In this section, we first discuss various aspects of the extended Kitaev model (including the three-spin interaction κ-term that breaks time-reversal symmetry) that are most relevant for this study.Second, we introduce the lattice Hamiltonian for the acoustic phonons on the honeycomb lattice and, third, perform the symmetry analysis of     (a) the magneto-elastic coupling.The explicit microscopic derivation of the coupling vertices in terms of Majorana fermions is later presented in Sec.IV.In Sec.III, we describe the dynamics of acoustic phonons both in the presence and in the absence of the time reversal symmetry using the long-wavelength effective action approach.
To obtain the hydrodynamic coefficients from the microscopic Hamiltonian, we use the diagrammatic techniques and compute the phonon polarization bubble.Details of the calculation of the bare Green's function and self-energy are presented in Sec.V.In Sec.VI, we relate the imaginary part of the diagonal components of the phonon polarization bubble to the attenuation coefficients.Phonon dynamics in the system with timereversal-symmetry breaking is discussed in Sec.VII.We first relate the off-diagonal component of the phonon polarization bubble to the Hall viscosity coefficient, and then show how it renormalizes the phonon spectrum.
Both the perturbative and non-perturbative corrections in terms of κ are obtained.A summary and a general discussion are given in Sec.VIII.Auxiliary information and technical details are provided in Appendices.

II. THE MODEL
We focus our discussion on the spin-phonon Hamiltonian The first term in Eq. ( 1) is the spin Hamiltonian given by where u α r,r+Mα = ib α r b α r+Mα , and r = r + M α and r = r+M α −M γ .Here we use Hs to label the Hamiltonian in the enlarged Hilbert space in terms of Majorana fermions.Importantly, the gauge dependent link variables u α r,r+Mα are mutually commuting constants of motion in Hs , and they give rise to physical static Z 2 flux degrees of freedom W p = Π r,r+Mα ∈p u α r,r+Mα = ±1 at the plaquette p.Each flux sector can be characterized by a static configuration of {u α r,r+Mα = ±1}, for which one obtains free-fermion Hamiltonian in terms of c r .Thus, the Majorana fermion Hamiltonian in the ground state zeroflux sector can be obtained by setting all u α r,r+Mα = 1.In momentum space, it reads where r∈B e ikr c r,B ) T are the complex fermions obtained from the Fourier transform of the Majorana fermions c A,r and c B,r on sublattices A and B, respectively, and N is the number of unit cells, τ x,y are the auxiliary Pauli matrices in sublattice space, ) from the nearest neighbor Majorana fermion hopping, and comes from the second neighbor Majorana hopping due to the time-reversal breaking part of the spin Hamiltonian, and 2 ).When κ = 0, the Majorana fermion spectrum k ∼ |f k | contains two Dirac points at the corners of the Brillouin zone ±K.
The second term in Eq. ( 1) is the bare Hamiltonian for the acoustic phonons on the honeycomb lattice, which contains the kinetic and elastic energy [51,52].The elastic part can be expressed in terms of the strain tensor ij = 1 2 (∂ i u j +∂ j u i ) and the elastic modulus tensor C ijlk , where u = {u x , u y } is the lattice displacement vector.From symmetry considerations, for a lattice with C 6v point group symmetry, there are only two independent non-zero elastic modulus tensor coefficients, C xxxx and C xxyy .The phonon Hamiltonian is written in terms of the bosonic operators that appear from quantizing the eigenmodes of the displacement vector u.We will postpone writing down the elastic energy until next section.
The third term in Eq. ( 1) denotes the magneto-elastic coupling that arises from the dependence of the Kitaev coupling J K on the relative positions between the spins.Assuming that J K only depends on the distance r between the atoms and that positions of the spins deviate only slightly from their equilibrium values, we expand the exchange couplings in powers of the displacement vectors.Keeping only the lowest order of this expansion, the spinphonon interaction part of the model Eq. ( 1) can then be written as where λ ∼ dJ K dr eq a characterizes the strength of the spin-phonon interaction and a is the lattice constant.
To find the complete set of magneto-elastic couplings, we use symmetry considerations and identify the linear combinations of the Kitaev interactions that transform as A sp 1 and E sp 2 irreducible representations (IRRs) of the C 6v point group.We find: Similarly, in the phonon sector, the strain combinations that transform as the IRRs A ph 1 and E ph 2 are ( xx + yy ) and { xx − yy , 2 xy }.Therefore, the spin-phonon couplings that are invariant under C 6v consist of two independent channels, one from A ph 1 ⊗ A sp 1 and another one from E ph 2 ⊗ E sp 2 , with the coupling constants λ A1 and λ E2 , which might be different but of a similar strength λ A1 , λ E2 ∼ dJ K dr eq a ∼ J K .Thus, the spin-phonon coupling Hamiltonian H c can be written as a sum of two independent contributions invariant under the C 6v symmetry, , where In Sec.IV, we will use the Majorana fermion representation of the spins to express the spin-phonon coupling H c in terms of the free Majorana fermion c r and lattice displacement field ũq,ν .

III. EFFECTIVE ACTION AND SPECTRUM FOR ACOUSTIC PHONON
In order to describe the dynamics of acoustic phonons both in the presence and in the absence of the time reversal symmetry breaking term κ, it is convenient to move away from the Hamiltonian formulation and employ instead the long-wavelength effective action S approach.Symmetries of the system, such as lattice symmetry and time-reversal symmetry, impose constraints on the number of non-zero independent coefficients in the elastic modulus tensor and viscosity tensor.In the following, only the non-dissipative terms in the effective action are considered.

A. Elastic medium with time-reversal symmetry
We start by considering acoustic phonons in a homogeneous space medium without decay and with timereversal symmetry.Their dynamics can be described by using the long-wavelength effective action in terms of the fields u, which describe the displacement of an atom from its original location.To lowest order, it reads [53] where ρ is the mass density of the lattice ions.Here the subscript in S (s) ph denotes that the elastic moduli C ijlk are symmetric under ij ↔ lk, and thus F (s) describes the symmetric part of the elastic energy.Due to the symmetry of the strain tensor ij under i ↔ j, the elastic modulus tensor C ijlk is also symmetric under i ↔ j and l ↔ k.Finally, by imposing the C 6v lattice symmetry we are left with two independent non-zero elastic modulus tensor coefficients, C xxxx and C xxyy .The elastic energy can then be written as: where In momentum and (Matsubara) frequency space, the action is given by We diagonalize Eq. ( 10) through u x,q u y,q = cos θ q − sin θ q sin θ q cos θ q ũq, ũq,⊥ , where ũq,ν denotes the magnitude of lattice displacement in the ν = , ⊥ channels which may be quantized as the acoustic phonons.We then find the longitudinal/transverse acoustic phonon spectrum and the polarization vectors (defined through u q = ν êq,ν ũq,ν ) to be equal to where q = q 2 x + q 2 y and θ q is defined in Fig. 2 (b).In the case of the bare lattice system, the phonon dispersions are isotropic, and the transverse/longitudinal modes are decoupled, which are generic properties for 2D medium with C 6 or higher symmetries.
Another property of elastic media is the viscosity tensor, which relates the stress tensor, obtained by differentiating the action with respect to the strain tensor, σ ij = − δS δ ij , to the strain rate, σ ij = −η ijkl ˙ kl .The viscosity tensor of non-interacting phonons is identically zero.However, it becomes finite when the phonon interactions that lead to phonon decay are included.In this work, we consider the non-zero viscosity tensor coefficients generated only through the spin-lattice coupling.The effects due to, e.g.phonon self-interactions, are ignored for simplicity.Moreover, generally, the correction from the phonon self-interactions is subleading in temperature [49].
As the viscosity tensor coefficient coming from the phonon interactions can be anisotropic in space, i.e. dependent on θ q , a straightforward symmetry analysis alone does not allow to determine the number of independent components of the symmetric viscosity tensor η (s) ijlk (by the same way as we just did for C ijlk ).As such, instead of finding the non-zero η (s) ijlk from a symmetry analysis, here we will calculate the decay rate directly from the phonon self-energy by integrating out the Majorana fermions.

B. Elastic medium with broken time-reversal symmetry
For an elastic medium with time-reversal symmetry breaking, e.g. when κ = 0 in Eq. ( 1), there is a direct contribution to the phonon effective action coming from the Hall viscosity term [48,49], where the viscosity tensor is anti-symmetric, i.e. η The lattice symmetries give additional constraints on the number of independent non-zero Hall viscosity tensor elements.In our case, there is only one independent component η (a) xxxy = η H .The reason is the following: the effective action must be invariant under all symmetry operations, so it must transform as the A 1 IRR.As the Hall viscosity tensor is anti-symmetric, we need to look for all the anti-symmetric A 1 IRRs formed by the tensor product of two phonon fields, i.e. those satisfying η ijlk = −η lkij .For the elastic medium with C 6v symmetry, two phonon fields in the E ph 2 irreducible representation (IRR) can form the tensor product given by 2 +..., where the superscript labels if the tensor is symmetric (s)/ anti-symmetric (a) under the exchange of the two phonon modes [54].With C 6v symmetry, there is no anti-symmetric A 1 IRR.On the other hand, the time-reversal symmetry breaking term κ in the spin Hamiltonian Eq. ( 2) lowers the symmetry from C 6v to C 6 .As a result, the A IRR of C 6 , following the compatibility relation [55].Consequently, we obtain one antisymmetric viscosity tensor element.
The effective action then becomes

S
(a) ph can also be expressed in terms of the transverse/longitudinal phonon eigenmodes ũq,ν (ν = , ⊥) as In general, the Hall viscosity coefficient η H can be found from the linear response theory through the Kubo formula, which gives [49] η µν (q) = lim for an action of the form where we remind that H c is the spin-phonon coupling Hamiltonian.
As we will show in Sec.VII A, the phonon polarization bubble Π µν ph (q, Ω) can be expressed in a similar fashion as Therefore, it is straightforward to find the relation between the response coefficient η µν (q) (and thus η H = 1 q 2 η ⊥ (q)) and Π µν ph (q, Ω).

IV. MICROSCOPIC DERIVATION OF THE EFFECTIVE LOW-ENERGY COUPLING HAMILTONIAN
In this section, we express the spin-lattice coupling in terms of the Majorana fermion-phonon (MFPh) coupling.To this end, we express the spin operators in Eq. ( 7) in terms of the Majorana fermions [3] and obtain the MFPh interaction in two symmetry channels: ), Since Majorana fermions satisfy {c r,α , c r ,β } = 2δ r,r δ α,β and {c k,α , c k ,β } = δ k,−k δ α,β , where α, β denote the sublattice indices A, B, the Fourier transformation of Majorana fermions is given by c r,A(B) = Thus, the MFPh coupling Hamiltonian in momentum space becomes with Here we used the fact that ij → i 2 (q i u j + q j u i ) and kept only leading in q terms.The auxiliary Pauli matrices τ x and τ y in sublattice space are again used in order to write the MFPh coupling in the matrix form.Note that we used a prime (double-prime) to denote the real (imaginary) part.We also defined In order to obtain the MFPh coupling vertices, we express the phonon modes in terms of the transverse/longitudinal eigenmodes defined in Eq.( 11).This gives where the MFPh vertices are To the leading order in momentum q and near the Dirac points ±K of the Majorana fermion spectrum, the MFPh vertices λν q,k are given by λ q,±K+k ≡ λ q,k =3 i q λ E2 (± sin 2θ q τ x − cos 2θ q τ y ), λ⊥ q,±K+k ≡ λ⊥ q,k =3 i q λ E2 (± cos 2θ q τ x + sin 2θ q τ y ), (28) where now k denotes the deviation from the Dirac point.Since the MFPh coupling in the A 1 channel appears at higher orders in k, whenever the computation is restricted to the low-energy Majorana fermions near the ±K points, the A 1 channel is ignored and only the dominant contribution in the E 2 channel is considered.

V. PHONON PROPAGATOR
Our next task is to integrate out the Majorana fermions and obtain an effective low-energy theory only for the phonons.

A. Free phonon propagator
The free phonon propagator in terms of lattice displacement field ũq,ν defined in Eq.( 11) is given by where the superscript (0) denotes the bare propagator, ν = , ⊥ labels the polarization.The second quantized form of ũq,ν (t) is given by [56] ũq,ν (t) = i 2ρ δ V Ω q 1/2 (a q e −iΩqt + a † −q e iΩqt ), (30) where δ V is the area enclosed in one unit cell and ρ is the mass density of the lattice ions.The time-ordered phonon propagator in the momentum and frequency space is then given by In the rest of the discussions, we set = 1.

B. Phonon polarization bubble
To compute the corrections to the effective phonon action due to the spin-lattice coupling, we calculate the phonon self-energy shown in Fig. 2 (b) that can be expressed as where µ, ν = , ⊥ and G(k, ω) denotes the Majorana fermions Green's function given by where ) is only non-zero when time-reversal symmetry is broken.We also use λq,k−q = λq,k +O(q 2 ), and drop the subleading O(q 2 ) terms.Tr[...] in Eq. ( 32) sums over momentum k, sublattice degrees of freedom (A and B) and integrates over frequency as dω 2π .Note that since the Majorana fermions are real, the Majorana fermionic excitation at k, ω and −k, −ω are physically the same, so G(k, ω) = −G(−k, −ω) T .We also note that since our computation is performed directly in terms of Majorana fermions, all k and ω modes should be summed over, which is different from the calculations that are performed in terms of the complex-fermion representation.
The renormalized phonon propagator is then given by the equation For calculations at finite temperature T , it is convenient to express the Majorana fermions Green's function and the phonon polarization bubble in the Matsubara frequency representation: and where Tr[...] now sums over momentum k, Matsubara frequencies iω m as T m and sublattice degrees of freedom.In general, the phonon polarization bubble Π ph contributes to the renormalization of the sound velocity [52], mixing of the transverse and longitudinal phonon modes [49], and the attenuation/absorption of sound waves [38,39].We compute these effects explicitly in Secs.VI and VII.As a final remark, we note that the corrections to the spin excitation dynamics from the spin-lattice coupling are relatively weak.The Majorana fermion self-energy [Fig. 2 (c)] can be expressed as where the bare phonon propagator in terms of the lattice displacement field ũq,ν is given by Eq.( 29).Here again we only keep terms linear in q and k in the MFPh coupling vertices, so λq,k−q λq,k .Tr[...] sums over momentum q, polarization ν and integrates over the phonon frequency as dΩ 2π .From Eq. ( 37), we find that the imaginary part of Σ ∼ T 2 , which is smaller than the typical fermion energy ∼ T .

VI. ANALYTICAL CALCULATION OF ATTENUATION COEFFICIENT
The quantitative description of the attenuation process can be obtained through the lossy acoustic wave function which decays with distance away from the driving source as u(x, t) = u 0 e −αs(q)x e i(Ωt−q•x) . ( where u(x, t) is the lattice displacement vector, u 0 = u(x, t = 0), Ω is the acoustic wave frequency and q is the propagation vector.The attenuation coefficient α s (q), defined as the inverse of the phonon mean free path, can be calculated from the imaginary part of the phonon selfenergy as where v s is the sound velocity.The derivation of this result is shown in Appendix B.
A. Sound attenuation coefficient at T > 0 The lowest-order Majorana fermion-phonon interaction leading to the phonon damping is shown diagrammatically in Fig. 2(b) with µ = ν.In this section, we will evaluate the sound attenuation coefficient α s (q) by relating it to the imaginary part of the diagonal components of the polarization bubble Im Π νν ph (q, Ω) defined in Eq.( 32).
We will consider the sound attenuation in the temperature range v F q < T < ∆ flux , where v F is the Majorana fermion velocity at the Dirac point, and assume v s < v F , the situation which can be potentially realized in Kitaev materials [8,46].Under these assumptions, there is a finite phase space for the scattering of phonons with thermally excited low-energy Majorana fermions [see Fig. 2  (a)].
We first note that from the kinematic constraints, Im Π νν ph (q, Ω) = 0 only at T > 0 when v s < v F .To compute Π νν ph (q, Ω) at finite temperature, it is convenient to first evaluate in the Matsubara frequency representation, i.e.Π νν ph (q, iΩ n ) in Eq. ( 36), and then perform the analytical continuation to the real frequency.Note that dynamical part only appears in the denominator of the Majorana fermion Green's function [see Eq. ( 35)].Thus, we can first evaluate the frequency summation [57]: where s, s correspond to the ± signs in the Majorana fermion Green's function G(k, iω m ) [see Eq. ( 35)].Performing the analytical continuation to the real frequency, we get Note that Im P s,s (k, q, Ω) = 0 only when there is a pole in the denominator of Eq. ( 41), i.e.Ω+s | F k−q |−s| F k | = 0, which for v s < v F requires the same choice of the sign in s and s , i.e. ss > 0.Moreover, when ss > 0, the numerator of Eq. ( 41) is finite only when In the temperature range of our interest, T < ∆ flux J K , the spectrum of the Majorana fermions can be linearized near the Dirac points ±K, and the MFPh vertices λν q,k are given by Eq. ( 28).In this approximation, the MFPh vertices are constant in terms of k and thus the momentum summation in Eq. ( 36) can be obtained by first replacing the true MFPh vertices with the Pauli matrices τ x,y .We then get where α, β = x, y.The explicit expressions for Πxx , Πyy and Πxy are derived in the Appendix C. The diagonal components of Im Π νν ph (q, Ω) can now be expressed in a straightforward way in terms of Παβ as Π ph (q, Ω) = 9λ 2 E2 q 2 [sin 2 2θ q Πxx + cos 2 2θ q Πyy − sin 2θ q cos 2θ q ( Πxy + Πyx )], Π ⊥⊥ ph (q, Ω) = 9λ 2 E2 q 2 [cos 2 2θ q Πxx + sin 2 2θ q Πyy + sin 2θ q cos 2θ q ( Πxy + Πyx )].(44) In the presence of time-reversal symmetry, i.e. when κ = 0, which corresponds to the z component of F k vanishing, we find that the imaginary part of the phonon polarization bubble for the longitudinal and transverse polarizations are (see Appendix C for derivation): In a system with time-reversal-symmetry breaking, from Eqs. ( 3) and ( 4), we note that the κ term does not excite Z 2 flux but opens a gap in the Majorana fermion bands with Chern number ±1 [3].The correction to the sound attenuation is perturbative in κ.Moreover, when the gap at ±K is small, i.e. κ < v F q, the correction is negligible.We do not consider this correction to the sound attenuation further.
To understand how the constraint v F q < T < ∆ flux is imposed in the computation, we first note that when v s < v F , the thermally excited Majorana fermions that can scatter off of the acoustic phonon should have minimum energy ∼ v F q/2.To optimize the phase space for phonon decay, we require v F q T .On the other hand, when T ∼ ∆ flux , the thermally excited flux degrees of freedom should be important, and the zero flux sector approximation breaks down.These processes remain to be considered in a future work.
From Eq. (B8), we find the sound attenuation coefficient to be equal to where 2 is determined by the elastic modulus tensor coefficient of the material, where C = C 1 + C 2 and C ⊥ = C 2 , and the upper/lower signs are for ν = and ν =⊥ polarizations, respectively.Importantly, (1 ∓ cos 6θ q ) is highly anisotropic and vanishes along certain directions, as shown in Fig. 1(b).
We note two important features of the attenuation coefficient.First, α s scales linearly with temperature, which counts the phase space of the Majorana fermion scattering.Second, while the phonon spectrum is isotropic at the leading q 2 order due to the six-fold rotation symmetry of the lattice, Im[Π νν ph (q, Ω)] can be anisotropic due to the coupling to the fermions with the Dirac spectrum.The angular modulation factor cos 6θ q may be understood from symmetry considerations as follows.The phonon polarization bubble [see Eq. ( 32)] is a convolution of two MFPh coupling vertices (each in the E ph 2 ⊗E sp 2 representation) and two Majorana fermion propagators [each G near the Dirac points in Eq. ( 33) contains in Eq. ( 32), which gives an angular modulation cos 6θ q .

VII. PHONON DYNAMICS WITH TIME-REVERSAL-SYMMETRY BREAKING
In this section, we will study the observable consequences of the Majorana fermion-phonon interaction in the presence of time reversal symmetry breaking, for example, due to an applied magnetic field.As was discussed above, the leading order perturbation from the magnetic field is the three-spin interaction κ term in the spin Hamiltonian, which does not change the sound attenuation coefficient qualitatively.On the other hand, the phonon system experiences Berry curvature induced by the κ term due to the spin-lattice coupling.In terms of the effective action, this is the Hall viscosity term given by Eq. ( 13).Here, we first compute the Hall viscosity coefficient by relating it to the off-diagonal component of the polarization bubble that is odd under time reversal, and then show how it renormalizes the phonon spectrum.

A. Hall viscosity coefficient with time-reversal symmetry breaking
We start by deriving the Hall viscosity coefficient η H which we introduced in Sec.III B. We first recall that it comes from the imaginary part of the off-diagonal component of the phonon polarization bubble, which is antisymmetric in exchanging the polarization indices µν.This process is non-dissipative as only the symmetric part of the viscosity coefficients contributes to dissipation [48][49][50]53].Consequently, Im Π ph (q, Ω) does not necessarily come from the pole in Π ph (q, Ω), which involves on-shell scattering and is generally dissipative.
In the following, we consider the contribution when on-shell processes are not involved, which is called the "intrinsic contribution" in the study of anomalous Hall effect [58].It turns out that the intrinsic contribution and the scattering contribution can be separated when v s < v F .Whereas the intrinsic contribution is non-zero already at T = 0, the scattering contribution requires T > 0 for the reason we discussed in Sec.VI A. We also note that as there is no kinetimatic constraint for the intrinsic contribution, both the low-energy and high-energy Majorana fermions contribute to the polarization bubble.
In the following, we restrict our analysis to T = 0. From Eqs. ( 32) and (33), we note that at T = 0 only the convolution of fermion propagators with ss < 0 are nonzero.Similarly to how the sound attenuation coefficient was obtained in the previous section, it is convenient to integrate over frequencies first.We obtain: Next, we compute the off-diagonal component of the polarization bubble.The MFPh coupling follows from Eqs. ( 26) and (27).To compute the polarization bubble, it is again convenient to first compute where τ α,β are Pauli matrices in sublattice space.Tr[...] sums over Matsubara frequency and sublattice degrees of freedom.We have The first (I) and the second (II) contributions have different physical meanings.The contribution (I) contains a term linear in Ω, and is odd under time-reversal, i.e.Ω → −Ω (in terms of both Matsubara frequency and time-ordered response in real frequency).Physically, this means that (I) contributes to the low-energy effective phonon action when time-reversal symmetry is broken.The contribution (II) renormalizes the real part of the phonon propagator, e.g. the sound velocity, which is a small effect and thus we do not consider further.
The above observations are manifest by computing Π ⊥ ph (q, iΩ) using the vertex functions λ (⊥)  q,k and summing k over the whole Brillouin zone.The MFPh coupling vertices originating from the time-reversal symmetry breaking spin Hamiltonian are ignored because they are smaller than those in Eq. ( 27) by a factor of κ/J K ∼ (h/J K ) 3 .Computing the contribution to Π ⊥ ph (q, iΩ) from part (I) in Eq. ( 50), which leads to the Hall viscosity term, we get where we have used Tr [τ α τ β τ γ ] = 2i αβγ , and defined λ (⊥) x,y as the τ x,y components of the MFPh coupling vertices λ (⊥)  q,k : Note that from Eq. ( 51) one clearly sees that the contribution to Π µν ph (q, iΩ) from the part (I) vanishes for µ = ν, and that Π ⊥ ph = −Π ⊥ ph .Furthermore, since ∆ k ∝ κ, the contribution from Eq. ( 51) only exists when time-reversal symmetry is broken (i.e. when κ = 0).
The main contribution to the momentum integration is near the Dirac points (±K).Linearizing the spectrum near the Dirac points gives 3 }, and Eq. ( 51) becomes We note that the result is "quantized" in the sense that it does not depend on microscopic details of the timereversal symmetry breaking perturbation, i.e. the value of κ, but only on its sign.In fact, one can show that it is essentially the same setting to compute Hall conductivity of anomalous quantum Hall effect when the MFPh vertex is constrained to near ±K [43].On the other hand, the analogy breaks down beyond the linear approximation of the Majorana fermion spectrum, since the high energy Majorana fermions also contribute to the Hall viscosity.
After integrating over the whole Brillouin zone, we find that the numerical prefactor depends on the magnitude of κ, i.e. on the magnitude of the Dirac mass.
In the following, we define the result found in Eq. ( 53) as Π ⊥ ph (q, iΩ) κ→0 , since in the limit |κ|/J K → 0 + , Eq. ( 53) is recovered.This is because the details of the Majorana fermion spectrum are not relevant in this limit, so the result should be the same as that obtained from the approximation of linearized spectrum.
The time-ordered phonon polarization bubble in real frequency can be obtained by analytic continuation.We have Here we need to stress that while both Eq. ( 54), leading to the Hall viscosity, and the imaginary part of the real time ordered phonon polarization bubble associated with the sound attenuation as we found in Eqs.(43) and (44) are linear in Ω, they are fundamentally different, i.e. the former is non-dissipative and the latter, dissipative.As was emphasized in our microscopic derivation, this difference comes from the distinct origins of the microscopic processes contributing to the imaginary part: one comes from the poles of the polarization bubble and another from the numerator ∼ Tr [τ x τ y τ z ].From a phenomenological perspective, the difference comes from the fact that the imaginary part of the time ordered polarization bubble associated with dissipation is even under Ω → −Ω due to causality, i.e. ∼ |Ω|.On the other hand, the imaginary part of the time ordered polarization bubble associated with the Hall viscosity is odd under Ω → −Ω, i.e. ∼ Ω.This is consistent with the causality constraint because the Hall viscosity term only appears in the off-diagonal component of the phonon polarization Π µν ph (q, Ω).The Hall viscosity coefficient η H in Eq. ( 16) can now be related to the phonon polarization bubble by expressing Eq. ( 32) in terms of where we recall that H c is the spin-lattice coupling Hamiltonian.Comparing with the expression in Eq. ( 16), we find the Hall viscosity coefficient to be equal to Finally, we also note that the Hall viscosity term may be the starting point to obtain the phonon thermal Hall conductivity, which is finite and of order ∼ T 3 [43,51].

B. Renormalization of the phonon spectrum
Next we discuss the renormalization of the phonon spectrum due to Π ⊥ ph (q, Ω) that mixes the longitudinal and transverse phonon modes.
For brevity in the expression, we write Eq. ( 54) as Im Π ⊥ ph (q, Ω) κ = q 2 Ω g(κ), where The renormalized phonon propagator can be expressed in matrix form as To the leading order in g(κ), the renormalized spectrum is given by Note that the correction from the q 3 term changes sign for the longitudinal and transverse phonon modes, so the spectrum bends upwards/downwards for the longitudinal/transverse mode [see Fig. 1 (c)].To determine the energy scale (relative to the Kitaev interaction J K ) when the deviation from the linear spectrum becomes prominent, the second term inside the brackets of the righthand side of Eq. ( 59) can be written as where the relation between the sound velocity, mass density and the components of the elastic modulus tensor was defined in Eq. ( 12), and It is straightforward to see that each bracket in Eq. (60) represents a dimensionless quantity.Thus, when Eq. ( 60) is comparable to 1, the cubic term in q becomes important.Therefore, when v ⊥ s < v s , the mixing shifts the transverse mode spectrum downwards and the longitudinal mode spectrum upwards at a characteristic scale

VIII. SUMMARY AND DISCUSSION
In this work, we studied the effects of the spin-lattice coupling on the acoustic phonon dynamics in the Kitaev spin liquid and obtained various observables, such as the sound attenuation (α s ), the Hall viscosity coefficient (η H ) and the renormalized 2D acoustic phonon spectrum [see Eq. ( 59)].We suggest that if measured, these observable effects can be used as potential probes of spin fractionalization in the Kitaev materials.
We demonstrated that these observables can be obtained from the phonon polarization bubble.To compute the bubble, a microscopic low-energy effective spin-lattice coupling was derived and a diagrammatic computation procedure was formulated in terms of the matter Majorana fermions and the acoustic phonons.We first showed that the sound attenuation comes from the phonon scattering off the Majorana fermions.Due to the temperature scaling of the scattering phase space, α s ∝ T , which distinguishes from the attenuation in other interaction channels, such as ∝ T 4 due to phonon self interactions.We also showed that α s is highly anisotropic with an angular modulation proportional to (1 ± cos 6θ q ), which should be attributed to the combined effects of the spinlattice coupling and the low-energy Dirac spectrum of the Majorana fermions.We then computed the Hall vis-cosity coefficient of the phonon effective action, which is non-zero when time-reversal symmetry is broken.Different from electron Hall fluids, the Hall viscosity coefficient is not quantized.We found two contributions: one that is non-perturbative and resembles the quantized Hall conductance, and another one that is perturbative in the time-reversal breaking term κ.The Hall viscosity coefficient may be probed indirectly from the phonon spectrum, where it is manifested as a deviation of the phonon dispersion from linear to non-linear near a characteristic momentum q c , which depends on the strength of the spin-lattice coupling [See Eq. ( 61)].Of course, such a deviation is generically expected for large enough momenta, where our long wavelength approximation breaks down.Our key point is that this new scale q c where the deviation from linear to non-linear takes place, is only present in the presence of an applied magnetic field.i.e. when time-reversal symmetry is broken.
As a final remark, we would like to emphasize that only in the temperature range T < ∆ flux , which we focus in this work, the thermal flux excitations can be ignored.On the other hand, the thermal flux excitations at higher temperatures may smear out some of the sharp features found in this work.Understanding the effects of spin-lattice coupling in the temperature range beyond the constraint set by ∆ flux would be highly desirable [40], and we leave it for future works.

FIG. 1 . 2 ) 3 2 , 1
FIG. 1.(a) The Kitaev model on the honeycomb lattice.The sites in the two sublattices of the honeycomb lattice (A and B) are marked by yellow and orange circles.n1 = ( √ 3 2 , 3 2 ) and n2 = (− √ 3 2 , 3 2 ) are two unit vectors, and we set the lattice constant a = 1.The three vectors Mx,y = (± √ 3 2 , 12 ), and Mz = (0, −1) connect nearest neighbors by x, y, and z bonds, respectively.The three sites r, r , r form a representative three-spin link r, r , r yx as defined below Eq. (2).(b) Angular dependence of the sound attenuation coefficient α s for the longitudinal mode.α ⊥ s for the transverse mode follows the same shape but is rotated by π/6.(c) Acoustic phonon spectrum Ω(q).The dashed lines depicts the bare spectrum with longitudinal and transverse modes.Solid lines depict the spectrum renormalized due to the mixing of transverse and longitudinal modes when time-reversal symmetry is broken.δc denotes the characteristic dimensionless scale beyond which the phonon spectrum is significantly bended [see Eq. (61)].
s h a 1 _ b a s e 6 4 = " I v N K Q e k r d r a V Q H U I k p u 7 I T C l S 1 w = " > A A A B 9 X i c b V B N S 8 N A E J 3 4 W e t X 1 a O X x S J 4 K o k I 9 l j w 4 s 0 K 9 g P a W D b b S b t 0 s 4 m 7 G 6 W E / g 8 v H h T x 6 n / x 5 r 9 x 2 + a g r Q y X p x 3 5 2 P e u u L k M 0 f w B 8 7 n D z 8 7 k k 0 = < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " I v N K Q e k r d r a V Q H U I k p u 7 I T C l S 1 w = " > A A A B 9 X i c b V B N S 8 N A E J 3 4 W e t X 1 a O X x S J 4 K o k I 9 l j w 4 s 0 K 9 g P a W D b b S b t 0 s 4 m 7 G 6 W E / g 8 v H h T x 6 n / x 5 r 9 x 2 + a g r Q s v L p H l e 8 d y K d 3 t R r l X z O A p w D C d w B h 5 c Q g 2 u o Q 4 N Y K D g G V 7 h z X l y X p x 3 5 2 P e u u L k M 0 f w B 8 7 n D z 8 7 k k 0 = < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " I v N K Q e k r d r a V Q H U I k p u 7 I T C l S 1 w = " > A A A B 9 X i c b V B N S 8 N A E J 3 4 W e t X 1 a O X x S J 4 K o k I 9 l j w 4 s 0 K 9 g P a W D b b S b t 0 s 4 m 7 G 6 W E / g 8 v H h T x 6 n / x 5 r 9 x 2 + a g r Q s v L p H l e 8 d y K d 3 t R r l X z O A p w D C d w B h 5 c Q g 2 u o Q 4 N Y K D g G V 7 h z X l y X p x 3 5 2 P e u u L k M 0 f w B 8 7 n D z 8 7 k k 0 = < / l a t e x i t > < l a t e x i t s h a 1 _ b a s e 6 4 = " I v N K Q e k r d r a V Q H U I k p u 7 I T C l S 1 w = " > A A A B 9 X i c b V B N S 8 N A E J 3 4 W e t X 1 a O X x S J 4 K o k I 9 l j w 4 s 0 K 9 g P a W D b b S b t 0 s 4 m 7 G 6 W E / g 8 v H h T x 6 n / x 5 r 9 x 2 + a g r Q s v L p H l e 8 d y K d 3 t R r l X z O A p w D C d w B h 5 c Q g 2 u o Q 4 N Y K D g G V 7 h z X l y X p x 3 5 2 P e u u L k M 0 f w B 8 7 n D z 8 7 k k 0 = < / l a t e x i t >

FIG. 2 .
FIG. 2. (a) On-shell process in which a phonon mode with Ωq = vs|q| (black dash arrow) excites a Majorana fermion from k to k + q such that ∆ = k+q − k = Ωq.Red (blue) shades indicate thermal excitation (depletion) probability of Majorana fermions near a Dirac cone at ±K.(b) Feynman diagram representing the phonon self-energy Eq. (36), which determines the sound attenuation and Hall viscosity for a given set of phonon polarization indices µν.(c) Feynman diagram representing the Majorana fermion self-energy Eq. (37) used to find the Majorana fermion lifetime due to spin-lattice coupling.

2 K.
)where J α denotes the nearest neighbor Kitaev interaction on the corresponding bond of type α = x, y, z, σ α r are the Pauli matrices, and M α labels the three inequivalent bonds on the honeycomb lattice [see Fig.1 (a)].In the following, we ignore the presence of weak non-Kitaev couplings and weak anisotropy of the Kitaev interactions in real materials, and assume J α = J K and the C 6v symmetry of the honeycomb lattice.The κ term breaks timereversal and mirror (σ v ) symmetry while preserving the exact solubility of the model.The three-spin link notation r, r , r αγ labels bonds rr , r r of type α, γ respectively on three adjacent sites r, r , r moving counterclockwisely [see Fig.1 (a)].α, β, γ are determined such that β = α, γ.It may be considered as the leading order perturbative effect of the magnetic field in the zero flux sector[3], and κ ∼ hxhyhz J Decomposing spin operators into two Majorana fermions σ α r = ic r b α r[3], the spin Hamiltonian can be written as lkij .In our case, the time-reversal symmetry breaking comes from the coupling of acoustic phonons with spin excitations, i.e. matter Majorana fermions with non-zero Chern number.After integrating out the fermion degrees of freedom, we are left with the phonon effective action S (a) ph that breaks time-reversal symmetry.

2
IRR of C 6v becomes the A (a) 1

2 N
k c k,A(B) e ik•r .In order to have simpler notations, in the following we denote c k,A ≡ a k and c k,B ≡ b k .