Topological magnon modes on honeycomb lattice with coupling textures

Topological magnon modes are expected to be useful for novel applications such as robust information propagation, since they are immune to backscattering and robust against disorder. Although there are several theoretical proposals for topological magnon modes and growing experimental efforts for realizing them by now, it is still desirable to add complementary insights on this important phenomenon. Here, we propose a new scheme to achieve topological magnon where only nearest-neighbour exchange couplings on honeycomb lattice are necessary. In both ferromagnets and antiferromagnets, tuning exchange couplings between and inside hexagonal unit cells induces a topological state accompanied by a band inversion between p-orbital and d-orbital like magnon modes. Topological magnon modes appear at the interface between a topological domain and a trivial domain with magnon currents, which counterpropagate depending on pseudospins originated from orbital angular momenta of magnon modes. This mimics the spin-momentum locking phenomenon in the quantum spin Hall effect.

n,k = i ψ ni b † i,k is the creation operator of n-th eigen magnon mode and ψ ni is the i-th element of ψ n . We notice that H F,k in Eq. (5) only differs from the electronic model on honeycomb lattice with hopping textures 35,36 , by the diagonal element E 0 , where t 0 and t 1 in electronic model are replaced by SJ 0 and SJ 1 respectively. Therefore, they share the same eigenstates with the shifted eigenvalues: where ε 0,n,k is the eigenvalue of H k in Eq. (9). One has ε F,n,k ∈ [0, 2E 0 ] , since ε 0,n,k ∈ [−E 0 , E 0 ]. The magnon frequency can be obtained by since ε F,n,k > 0.
The frequency band structures of magnon modes are shown in Fig. 2. For J 0 = J 1 , the frequency band structure exhibits Dirac cones same as that of the electronic structure of graphene 30 . These Dirac cones are folded to the Ŵ point as shown in Fig. 2b, since we are using the six-site hexagonal unit cell instead of the conventional two-site rhombic unit cell. For J 0 = J 1 , a frequency gap opens at the Ŵ point as displayed in Fig. 2a,c. As clarified in Ref. 32 , the band structures for J 0 > J 1 (Fig. 2a) and for J 0 < J 1 (Fig. 2c) are distinct in topology.
Because Hamiltonian (1) preserves the inversion symmetry, eigen wavefunctions of Hamiltonian (5) can be indexed by parity eigenvalues, even and odd, at the Ŵ point and M point. In the case J 0 > J 1 (see Fig. 2a), the two degenerate wavefunctions below the gap are |p x � and |p y � with odd parity at the Ŵ point, whereas the two degenerate wavefunctions above the gap are |d x 2 −y 2 � and |d 2xy � with even parity. The p and d magnon modes are shown in Fig. 2d. The parity of each band at the M point is the same as that at the Ŵ point, denoting a topologically trivial  38,39 . In the case J 0 < J 1 , the frequencies of p and d modes are inverted at the Ŵ point as shown in Fig. 2c. Now, at the Ŵ point, all three eigenstates below the gap are of even parity, while at the M point, there are two eigenstates with odd parity and one eigenstate with even parity. The unequal numbers of eigenstates with even parity indicate a topological state 38,39 . The band structure of topological state in Fig. 2c cannot be continually transformed from that of the topologically trivial state in Fig. 2a without closing the gap (see Fig. 2b). The double degeneracy of p (d) modes at the Ŵ point are protected by the C 6v symmetry, where the wavefunctions of p (d) modes correspond to the basis functions of two-dimensional irreducible representation E 1 ( E 2 ). In order to better understand the symmetry at the Ŵ point, we introduce a pseudo time-reversal operator T = iσ y K , where σ y is Pauli matrix referring to the basis functions of E 1 and E 2 , and K denotes the complex conjugate operator 32,35 . This pseudo time-reversal operator originates from the C 6v symmetry and the time-reversal symmetry of the system, which produces the Kramers doubling in the present bosonic system. Kramers pairs associated with p or d orbitals are constructed as |p ± � = 1 √ 2 (|p x � ± i|p y �) and |d ± � = 1 √ 2 (|d x 2 −y 2 � ± i|d 2xy �) respectively, which carry specific pseudospins defined on the "artificial molecule". In the present system, the magnon current I ij = �Î ij � is given by the current operator For magnon modes with pseudospin up/down, magnon currents circulate counterclockwise/clockwise, as shown in Fig. 2e,f. (|p x � ± i|p y �) and |d ± � = 1 √ 2 (|d x 2 −y 2 � ± i|d 2xy �). www.nature.com/scientificreports/ Around the Ŵ point, magnon modes are predominantly occupied by the p orbitals and d orbitals. We can rewrite Hamiltonian (5) in the basis of [|p + �, |d + �, |p − �, |d − �] up to the lowest order of momentum k as with where a 0 is the length of unit vector (see Fig.1). This effective Hamiltonian is similar to the Bernevig-Hughes-Zhang model for quantum spin Hall effect and that found for topological photonic crystals 3,32 . In Hamiltonian (17), E 0 simply makes a constant shift in the eigenvalues, and is the effective Dirac mass, which is positive in a topologically trivial state. In the case J 0 < J 1 , M becomes negative and a band inversion takes place, which turns the system into a topological state with a topological gap 2(J 1 − J 0 ) , as discussed above and shown in Fig. 2c. It is noticed that the topological magnon modes in the present system rely on crystalline symmetry and have a weaker topology comparing to those systems with the DM interaction. Nevertheless, we expect that stable unidirectional interface states between topological and trivial domains can be observed experimentally, similarly to phononic and photonic systems [40][41][42][43] .
In order to investigate topological interface magnon modes, we consider a heterostructure as illustrated in Fig. 3a, which is uniform and infinitely long in the x direction. In the y direction, the heterostructure contains 60 unit cells in both trivial and topological domains periodically, for the simplicity of calculation. Exchange couplings inside trivial and topological domains are of the same values used in Fig. 2a,c respectively, whereas exchange couplings between the trivial and topological domains, which are denoted by gray lines in Fig. 3a, are given by the geometric mean of value of J 1 used in Fig. 2a,c. For each interface, two interface dispersions appear inside the bulk gap in the frequency band structure shown in Fig. 3b, which is obtained numerically using the large unit cell for the heterostructure (Fig. 3a). The Hamiltonian of heterostructure takes the same form as Hamiltonian (5) except for that it is for a supercell (see Fig. 3a) and the momentum is along the x direction. The Dirac mass takes opposite signs in the topological and trivial domains, leading to Jackiw-Rebbi soliton solutions localized at the interface, whose amplitudes decay exponentially into the bulks 44 . The interface magnon modes at the two momenta marked by ① and ② in Fig. 3b are shown in Fig. 4a,b respectively. The topological magnon modes correspond to precessions of spins around the z axis, with amplitudes decaying into the bulks exponentially. These magnon modes can be excited by applying an external oscillating magnetic field. In this way, the frequency of excited magnon mode is controllable. In Fig. 4c,d, magnon currents of topological interface magnon modes with counterclockwise circulation and clockwise circulation in unit cells correspond to the upand down-pseudospin respectively. The net current for up-/down-pseudospin flows to the positive/negative x direction, manifesting the pseudospin-momentum locking in the present topological magnon modes.
Topological magnon modes in antiferromagnets. Next, we explore topological magnon modes in antiferromagnets on honeycomb lattice. Antiferromagnets exhibit degenerate dispersions guaranteed by the combination of time-reversal and inversion symmetries, where no Dirac cone exists 33 . Therefore, it is not straightforward to realize topological magnon modes using the procedure formulated for ferromagnets directly. In order to overcome this difficulty, we introduce an external magnetic field to lift the double degeneracy. The Heisenberg model for antiferromagnets is where exchange couplings J ij > 0 , the external magnetic field B ext > 0 , and spin lengths S i = S , and for simplicity, the external magnetic field is perpendicular to the two-dimensional spin plane. Without single ion anisotropy considered here, spins are automatically aligned in the plane perpendicular to the direction of external magnetic field and canted as shown schematically in Fig. 5a. The canting angle θ of an antiferromagnet under external magnetic field is determined by the competition between exchange couplings and the Zeeman energy. We consider uniform exchange couplings J ij = J as shown in Fig. 5a. The energy of a classical canted antiferromagnetic state with N rhombic unit cells is which is minimized when spins are canted at an angle sinθ = B ext /(6JS).
In order to apply the Holstein-Primakoff transformation for antiferromagnets, we need a rotating frame where the z coordinate axis is rotated to the local spin direction at each site (see Fig. 5b). The transformation between the rotating frame and the laboratory frame is where S x , S y , S z are spin components of the laboratory frame, S ′x , S ′y , S ′z are those of the rotating frame, φ is the angle between the projection of a spin in the xy plane and the coordinate axis x and θ is the canting angle of spin. Note that projections of spins on different sublattices point to opposite directions on the xy plane, namely with δ i being the n.n. vectors shown in Fig. 5a. Note that the terms linear in b i or b † i cancel out for sinθ = B ext /(6JS) , which minimizes Eq. (20). Because the classical Néel order is not the ground state for quantum antiferromagnets, bb and b † b † terms that do not conserve the total number of magnon appear, as such the basis vector k including both creation and annihilation operators is chosen.
Dispersions of magnon modes can be derived by applying the bosonic Bogoliubov transformation 45,46 T is a new set of creation and annihilation operators for magnons and E k is a diagonalized matrix. Magnons both before and after the Bogoliubov transformation should obey the same commutation relation for bosons in Eq. (12), demanding that the matrix T k satisfies Therefore, one cannot calculate E k from the characteristic equation of matrix H CAF,k directly, since columns of T k are not eigenvectors of H CAF,k .
Multiplying Î (T † k ) −1 to both sides of T † k H CAF,k T k = E k from left, we obtain the relation namely columns of T k are eigenvectors of Î H CAF,k . The diagonalized matrix Î E k can be derived from the characteristic equation of matrix Î H CAF,k , det|ÎH CAF,k −ÎE k | = 0 , which is the conventional way to obtain Î E k and T k . Because k and Ŵ k include both creation and annihilation operators of magnons, dispersions in Eq. (25) are redundant, which should be eliminated by symmetries. First, we define a particle-hole operator where I 2 is a 2 × 2 identity matrix, K is the complex conjugate operator and P is the space-inversion operator. For Î H CAF,k , we can directly check from Eqs. (24) and (27) that Î H CAF,k C = −CÎH CAF,k , leading to Î E k = diag(ε 1,k , ε 2,k , −ε 1,k , −ε 2,k ) . With the same method, we can also check that PÎH CAF,k =ÎH CAF,−k P , which leads to ε n,k = ε n,−k ( n = 1, 2 ). In this way, Hamiltonian (23) is diagonalized into [45][46][47] with where (u A,n,k , u B,n,k , v * A,n,−k , v * B,n,−k ) is the nth column of T k . Then the dispersions of magnon modes for the canted antiferromagnet shown in Fig. 5a can be derived from the above procedure explicitly as 33 It is clear that without an external magnetic field (θ = 0) , the dispersion is doubly degenerate in the whole Brillouin zone (see Fig. 5c) where no Dirac cone exists, unlike ferromagnets (see Fig. 2b). A finite external (28)  www.nature.com/scientificreports/ magnetic field lifts the degeneracy, except for the K and K ′ points where f K = f K ′ = 0 and ε i,K = ε i,K ′ = 3SJ , as shown in Fig. 5d. By choosing a hexagonal unit cell for canted antiferromagnets with a coupling texture as shown in Fig. 6a, where the canting angle is determined by sinθ = B ext /2S(2J 0 + J 1 ) now, Hamiltonian (24) becomes where E 0 = S(2J 0 + J 1 ) , I 6 is a 6 × 6 identity matrix and H F,k is Hamiltonian (5) for ferromagnets.
Eigenvalues of canted antiferromagnets can be derived in terms of eigenvalues of ferromagnets with the same coupling texture. Here we consider the particle solution of Eq. (26) which has positive energy ε n,k with n = 1, 2, ..., 6. All four blocks of H CAF,k in Eq. (31) are digonalized in the basis of eigenstates of H F,k , so that for an eigenstate ψ of H F,k with eigenvalue ε F,n,k , the state (ψ, αψ) T is an eigenstate of Î H CAF,k in Eq. (26), provided leading to For sin 2 θ > 1/3 , ε n,k is monotonic with ε F,n,k . The p-d band inversion is achieved in canted antiferromagnets by tuning coupling texture the same as in the ferromagnetic case. For sin 2 θ = 0.4 , the full frequency band structures of canted antiferromagnets with the coupling textures are derived as shown in Fig. 6b-d, and p and d magnon modes in canted antiferromagnets are shown in Fig. 6e. For sin 2 θ < 1/3 , the present scheme does not apply straightforwardly.
In order to see the topological interface magnon modes, we consider a heterostructure shown in Fig. 7a similar to Fig. 3a for the ferromagnet. Topological interface dispersions appear inside the band gap of the frequency band structure as shown in Fig. 7b. Magnon modes at the two momenta marked by ① and ② in Fig. 7b are depicted in Fig. 8a,b respectively. We can calculate the magnon current for particle and hole parts of canted antiferromagnet using the same way of ferromagnetic case. The total magnon current is the magnon current of the particle part minus that of the hole part. Similar to topological interface magnon modes of the ferromagnetic case, in Fig. 8c,d magnon current circulates counterclockwise/clockwise in unit cells dominated by the up-/down-pseudospin, which also governs the direction of net magnon currents, demonstrating the pseudospin-momentum locking  Fig. 5a, where the rotating frame at A and B sublattice are different by a π rotation of φ , magnon modes look exactly the same way as Fig. 2d. For simplicity, φ = ±π/2 is taken. www.nature.com/scientificreports/ phenomenon in these topological interface magnon modes. Because the external magnetic field induces a ferromagnetic component Ssinθ on each site, the magnitude of magnon current is proportional to sin 2 θ.

Discussion
We propose a method to achieve topological magnon modes in magnetic systems on honeycomb lattice, including both ferromagnet and antiferromagnet. The frequency band structures are gapless for uniform nearest-neighbor exchange couplings. In ferromagnets, a topological frequency gap opens when exchange couplings inside the hexagonal unit cells are smaller than exchange couplings between unit cells, associated with a p-d band inversion at the Ŵ point. In antiferromagnets, the degeneracy in the frequency band structure due to the combination of time-reversal symmetry and inversion symmetry has to be lifted by applying an external magnetic field. The resulting canted antiferromagnets become topological upon tuning exchange couplings inside and between hexagonal unit cells same as in ferromagnets.
In the present work, a hexagonal unit cell with six sites instead of a rhombic cell with two sites is used. It is crucially important in our proposal since the present crystalline topology is based on the C 6v symmetry which should be preserved even when the coupling texture is introduced. On the other hand, the conditions 2J 0 + J 1 = 3J and |J 0 − J 1 | = 0.3J are considered for simplicity for which the band-gap center and size are the same for topological and trivial magnon states. Topological interfacial magnon modes can be found so long as J 0 > J 1 and J 0 < J 1 are satisfied in trivial and topological domains respectively and their bulk gaps overlap partially. In the present work, only short-range exchange interactions are considered in the Hamiltonian, whereas long-range dipolar effects are neglected.
Magnon currents of topological magnon modes propagate along the interface between a trivial domain and a topological domain in opposite directions governed by pseudospins, the circulation direction of magnon current in unit cells, manifesting the pseudospin-momentum locking phenomenon. Interesting phenomena may be caused by higher order magnon terms such as the 4-magnon scattering on the thin layers 48 , but they are beyond the scope of our present study which focusses on the linear part of the magnon Hamiltonian.  Fig. 6b,d respectively, which are denoted in the same way as in Fig. 3a, and the couplings between the trivial and topological domains (denoted by gray lines) are given by geometric mean of J 1 . (b) Frequency band structure of magnon modes for the heterostructure, where topological interface dispersions (red lines) appear in the bulk band gap. Wavefunctions and magnon currents of magnon modes at the two opposite momenta denoted by 1 and 2 will be depicted in Fig. 8 www.nature.com/scientificreports/ In the future, candidate materials should be found to realize the scheme proposed in the present work. Besides specific materials, one may also observe topological magnon modes in artificial systems which can be achieved by depositing magnetic atoms on a metallic substrate using the STM technique or trapping magnetic atoms in an optical lattice using laser beams. In these artificial systems, the tuning of exchange couplings proposed in this work may be realized by manipulating the distance of neighbor magnetic atoms.

Methods
The effective Hamiltonians of magnon modes in ferromagnet and antiferromagnet on honeycomb lattice are obtained by the Holstein-Primakoff approach. Besides, a local rotating frame and the Bogoliubov transformation are adopted for canted antiferromagnet. Frequency band structures, wavefunctions of eigen modes and topological interface modes in the present work, such as those in Figs. 2a-d and 3b, are obtained by diagonalizing Hamiltonians of bulk and heterostructure, respectively, using MATLAB.