Three-dimensional topological magnon systems

TL;DR

This paper introduces three-dimensional topological magnon models with pseudo-time-reversal symmetry, defining $Z_2 invariants, and demonstrates phases with surface Dirac cones, suggesting potential thermal Hall effects in magnetic systems.

Contribution

It develops a bosonic topological insulator model with $Z_2 invariants and explores its phases and surface states, extending topological insulator concepts to magnonic systems.

Findings

Identification of three topological phases in magnon systems

Presence of surface Dirac cones in non-trivial phases

Potential for thermal Hall effect in surface magnons

Abstract

We propose a class of models for a magnonic analog of topological insulators in three dimensions. The models have pseudo-time-reversal symmetry which ensures the existence of bosonic Kramers pairs. We define a set of $\mathbb{Z}_2$ topological invariants that characterizes different topological phases and determines the presence or absence of surface Dirac cones. This is demonstrated by considering a bosonic counterpart of the Fu-Kane-Mele model on a diamond lattice. The model is found to exhibit three distinct phases analogous to strong topological, weak topological, and trivial insulator phases of the original fermionic model. We also discuss a possible realization of the thermal Hall effect of surface magnons in the presence of a magnetic field in proximity to a normal ferromagnet.

Tables

Table 1. Table 1 : The ℤ 2 subscript ℤ 2 \mathbb{Z}_{2} topological invariants of the model. The parameters are chosen to be J 0 = 1.4 , J 1 = J 2 = J 3 = J ′ = 1.0 , J − = D = Γ = 0.3 , κ = 1.5 formulae-sequence formulae-sequence subscript 𝐽 0 1.4 subscript 𝐽 1 subscript 𝐽 2 subscript 𝐽 3 superscript 𝐽 ′ 1.0 subscript 𝐽 𝐷 Γ 0.3 𝜅 1.5 J_{0}=1.4,J_{1}=J_{2}=J_{3}=J^{\prime}=1.0,J_{-}=D=\Gamma=0.3,\kappa=1.5 . The index n = 1 , 2 𝑛 1 2 n=1,2 denotes the top and bottom band in the particle space, respectively.

$n$	1	2
${\nu }_{x,0}^n$	0	0
${\nu }_{x,\pi }^n$	1	1
${\nu }_{y,0}^n$	0	0
${\nu }_{y,\pi }^n$	1	1
${\nu }_{z,0}^n$	0	0
${\nu }_{z,\pi }^n$	1	1
$({\nu }_0^n;{\nu }_x^n,{\nu }_y^n,{\nu }_z^n)$	(1;1,1,1)	(1;1,1,1)

Table 2. Table 2 : The topological invariants of the perovskite lattice system. Parameters are chosen to be J ​ S = 0.03 , D ​ S = 0.4 , J + ​ S = 0.1 , J − ​ S = 0.1 , J ′ ​ S = 1.0 , κ A ​ S = 1.35 , κ B ​ S = 0.5 , κ C ​ S = 1.0 , s 1 = − ( 2 + 2 ) / 4 , s 2 = ( 2 − 2 ) / 4 , s 3 = 1 / 2 formulae-sequence 𝐽 𝑆 0.03 formulae-sequence 𝐷 𝑆 0.4 formulae-sequence subscript 𝐽 𝑆 0.1 formulae-sequence subscript 𝐽 𝑆 0.1 formulae-sequence superscript 𝐽 ′ 𝑆 1.0 formulae-sequence subscript 𝜅 𝐴 𝑆 1.35 formulae-sequence subscript 𝜅 𝐵 𝑆 0.5 formulae-sequence subscript 𝜅 𝐶 𝑆 1.0 formulae-sequence subscript 𝑠 1 2 2 4 formulae-sequence subscript 𝑠 2 2 2 4 subscript 𝑠 3 1 2 JS=0.03,DS=0.4,J_{+}S=0.1,J_{-}S=0.1,J^{\prime}S=1.0,\kappa_{A}S=1.35,\kappa_{B}S=0.5,\kappa_{C}S=1.0,s_{1}=-(2+\sqrt{2})/4,s_{2}=(2-\sqrt{2})/4,s_{3}=1/2 . The index n = 1 , 2 , 3 𝑛 1 2 3 n=1,2,3 denotes the band with n 𝑛 n th highest energy in the particle space.

$n$	${\nu }_{x,0}^n$	${\nu }_{x,\pi }^n$	${\nu }_{y,0}^n$	${\nu }_{y,\pi }^n$	${\nu }_{z,0}^n$	${\nu }_{z,\pi }^n$	$({\nu }_0^n;{\nu }_x^n,{\nu }_y^n,{\nu }_z^n)$
1	1	0	0	1	0	1	(1;0,1,1)
2	1	0	0	1	0	1	(1;0,1,1)
3	0	0	0	0	0	0	(0;0,0,0)

Table 3. Table 3 : The topological invariants of the pyrochlore lattice system. Parameters are chosen to be J ​ S = 0.1 , D ~ ​ S = 0.4 , J + ​ S = 0.8 , J − ​ S = 0.6 , J ′ ​ S = 1.0 , κ A ​ S = 3.0 , κ B , C , D ​ S = 4.5 formulae-sequence 𝐽 𝑆 0.1 formulae-sequence ~ 𝐷 𝑆 0.4 formulae-sequence subscript 𝐽 𝑆 0.8 formulae-sequence subscript 𝐽 𝑆 0.6 formulae-sequence superscript 𝐽 ′ 𝑆 1.0 formulae-sequence subscript 𝜅 𝐴 𝑆 3.0 subscript 𝜅 𝐵 𝐶 𝐷 𝑆 4.5 JS=0.1,\tilde{D}S=0.4,J_{+}S=0.8,J_{-}S=0.6,J^{\prime}S=1.0,\kappa_{A}S=3.0,\kappa_{B,C,D}S=4.5 . The index n = 1 , 2 , 3 , 4 𝑛 1 2 3 4 n=1,2,3,4 denotes the band with the n 𝑛 n th highest energy in the particle space.

$n$	${\nu }_{x,0}^n$	${\nu }_{x,\pi }^n$	${\nu }_{y,0}^n$	${\nu }_{y,\pi }^n$	${\nu }_{z,0}^n$	${\nu }_{z,\pi }^n$	$({\nu }_0^n;{\nu }_x^n,{\nu }_y^n,{\nu }_z^n)$
1	1	1	1	1	1	1	(0;1,1,1)
2	0	0	0	0	0	0	(0;0,0,0)
3	1	0	1	0	1	0	(1;0,0,0)
4	0	1	0	1	0	1	(1;1,1,1)

Equations

$$\displaystyle\mathcal{H}=\frac{1}{2}\sum_{\bm{k}}[\bm{\beta}^{\dagger}(\bm{k})\bm{\beta}(-\bm{k})]H(\bm{k})\left[\begin{array}[]{cc}\bm{\beta}(\bm{k})\\ \bm{\beta}^{\dagger}(-\bm{k})\\ \end{array}\right].$$

$$⟨⟨\mathbf{\varphi },\mathbf{\psi }⟩⟩={\mathbf{\varphi }}^{†}{\Sigma }_z\mathbf{\psi },$$

$${\Sigma }_{z}H(-\mathbf{k}){\Theta }^{\prime }-{\Theta }^{\prime }{\Sigma }_{z}H(\mathbf{k})=0.$$

$${\nu }_{i,0}^{n\sigma }:=\frac{1}{2\pi }\left[{\oint }_{\partial {\mathrm{EBZ}}_{i,0}}d\mathbf{k}\cdot {\left[{\mathbf{A}}_{n\sigma }(\mathbf{k})\right]}_{k_i=0}-{\int }_{{\mathrm{EBZ}}_{i,0}}d{k}_{j}d{k}_k{\left[{\Omega }_{n\sigma }^i(\mathbf{k})\right]}_{k_i=0}\right]\mathrm{mod}2,$$

$${\nu }_{i,\pi }^{n\sigma }:=\frac{1}{2\pi }\left[{\oint }_{\partial {\mathrm{EBZ}}_{i,\pi }}d\mathbf{k}\cdot {\left[{\mathbf{A}}_{n\sigma }(\mathbf{k})\right]}_{k_i=\pi }-{\int }_{{\mathrm{EBZ}}_{i,\pi }}d{k}_{j}d{k}_k{\left[{\Omega }_{n\sigma }^i(\mathbf{k})\right]}_{k_i=\pi }\right]\mathrm{mod}2,$$

$${\mathbf{A}}_{n\sigma }(\mathbf{k})=\sum _l{\mathbf{A}}_{nl\sigma }(\mathbf{k}),$$

$${\mathbf{\Omega }}_{n\sigma }(\mathbf{k})=\sum _l{\mathbf{\Omega }}_{nl\sigma }(\mathbf{k}),$$

$${\mathbf{A}}_{nl\sigma }(\mathbf{k})=\mathrm{i}⟨⟨{\mathbf{\Psi }}_{nl\sigma }(\mathbf{k}),{\nabla }_{\mathbf{k}}{\mathbf{\Psi }}_{nl\sigma }(\mathbf{k})⟩⟩,$$

$${\mathbf{\Omega }}_{nl\sigma }(\mathbf{k})={\nabla }_{\mathbf{k}}\times {\mathbf{A}}_{nl\sigma }(\mathbf{k}).$$

$${\nu }_{i,0}:=\sum _{n,(E_n(\mathbf{k})\le {ϵ}_0)}{\nu }_{i,0}^n\mathrm{mod}2,$$

$${\nu }_{i,\pi }:=\sum _{n,(E_n(\mathbf{k})\le {ϵ}_0)}{\nu }_{i,\pi }^n\mathrm{mod}2,$$

$${\nu }_{x,0}+{\nu }_{x,\pi }={\nu }_{y,0}+{\nu }_{y,\pi }={\nu }_{z,0}+{\nu }_{z,\pi }\mathrm{mod}2.$$

$${\nu }_0={\nu }_{x,0}+{\nu }_{x,\pi }\mathrm{mod}2,$$

$${\nu }_i={\nu }_{i,\pi }(i=x,y,z).$$

$${\mathbf{\beta }}^{†}(\mathbf{k})=[{\mathbf{b}}_{\uparrow }^{†}(\mathbf{k}),{\mathbf{b}}_{\downarrow }^{†}(\mathbf{k})].$$

$${\mathbf{b}}_{\uparrow }^{†}(\mathbf{k})=[b_{\uparrow ,1}^{†}(\mathbf{k}),\cdots ,b_{\uparrow ,N}^{†}(\mathbf{k})],$$

$${\mathbf{b}}_{\downarrow }^{†}(\mathbf{k})=[b_{\downarrow ,1}^{†}(\mathbf{k}),\cdots ,b_{\downarrow ,N}^{†}(\mathbf{k})],$$

$${\Theta }^{\prime }=({\sigma }_z\otimes \mathrm{i}{\sigma }_y\otimes 1_N)K.$$

$$H=H_{\mathrm{DM}}+H_{J^{\prime }}+H_J+H_{\mathrm{XY}}+H_{\mathrm{\Gamma }}+H_{\kappa }.$$

$$H_J\leftrightarrow t\sum _{⟨i,j⟩}{c}_i^{†}{c}_j,$$

$$H_{\mathrm{DM}}\leftrightarrow i(8{\lambda }_{SO}/a^2)\sum _{⟨⟨i,j⟩⟩}{c}_i^{†}{s}_z({\mathbf{d}}_{ij}^1\times {\mathbf{d}}_{ij}^2{)}_z{c}_j,$$

$$H_{\mathrm{\Gamma }}\leftrightarrow i(8{\lambda }_{SO}/a^2)\sum _{⟨⟨i,j⟩⟩}{c}_i^{†}{s}_x({\mathbf{d}}_{ij}^1\times {\mathbf{d}}_{ij}^2{)}_x{c}_j,$$

$$H_{\mathrm{XY}}\leftrightarrow i(8{\lambda }_{SO}/a^2)\sum _{⟨⟨i,j⟩⟩}{c}_i^{†}{s}_y({\mathbf{d}}_{ij}^1\times {\mathbf{d}}_{ij}^2{)}_y{c}_j.$$

$$\displaystyle H(\bm{k})=\left(\begin{array}[]{cc}h(\bm{k})&\Delta^{\dagger}(\bm{k})\\ \Delta(\bm{k})&h^{*}(-\bm{k})\\ \end{array}\right),$$

$$h(\mathbf{k})=d_0{1}_4+\sum _{a=1}^5{d}_a(\mathbf{k}){\Gamma }^a,$$

$$\Delta (\mathbf{k})=J^{\prime }{\sigma }_x\otimes 1_2,$$

$${\Gamma }^{(1,2,3,4,5)}=(1_2\otimes {\sigma }_x,1_2\otimes {\sigma }_y,{\sigma }_x\otimes {\sigma }_z,{\sigma }_y\otimes {\sigma }_z,{\sigma }_z\otimes {\sigma }_z).$$

$$d_0=J_0+J_1+J_2+J_3+J^{\prime }+2\kappa ,$$

$$d_1(\mathbf{k})=-J_0-J_1\cos (k_1)-J_2\cos (k_2)-J_3\cos (k_3),$$

$$d_2(\mathbf{k})=-J_1\sin (k_1)-J_2\sin (k_2)-J_3\sin (k_3),$$

$$d_3(\mathbf{k})=-\sqrt{2}\Gamma [\sin (k_2)-\sin (k_3)-\sin (k_{21})+\sin (k_{31})],$$

$$d_4(\mathbf{k})=-\sqrt{2}{J}_{-}[\sin (k_3)-\sin (k_1)-\sin (k_{32})+\sin (k_{12})],$$

$$d_5(\mathbf{k})=-\sqrt{2}D[\sin (k_1)-\sin (k_2)-\sin (k_{13})+\sin (k_{23})],$$

$$\Delta E=\sqrt{(d_0+|d_1(\mathbf{\Gamma })|{)}^2-J^{\prime 2}}-\sqrt{(d_0-|d_1(\mathbf{\Gamma })|{)}^2-J^{\prime 2}},$$

$$|\mathbf{\Psi }(\mathbf{k})⟩$$

$$:=\left[|{\mathbf{\Psi }}_{n,1,+}(\mathbf{k})⟩,|{\mathbf{\Psi }}_{n,2,+}(\mathbf{k})⟩,|{\mathbf{\Psi }}_{n+1,1,+}(\mathbf{k})⟩,|{\mathbf{\Psi }}_{n+1,2,+}(\mathbf{k})⟩\right].$$

$$H_n(\mathbf{k})=⟨\mathbf{\Psi }(\mathbf{k})|H(\mathbf{k})|\mathbf{\Psi }(\mathbf{k})⟩.$$

$$H_z(-\mathrm{i}{\partial }_z)\mathbf{\psi }(z)=E_0\mathbf{\psi }(z).$$

Full text

I INTRODUCTION

In recent years, there has been growing interest in studying topological phases of matter Schnyder08 ; Kitaev09 ; Ryu10 . In particular, the discovery of quantum spin Hall insulators Kane05a ; Kane05b attracted a great deal of attention. The topological phases of quantum spin Hall insulators are characterized by the ${\mathbb{Z}}_{2}$ topological invariant Fu06 ; Fukui07 ; Fukui08 ; Qi08 ; Wang10 ; Loring10 ; Fulga12 ; Sbierski14 ; Loring15 ; Katsura16 ; Akagi17 ; Katsura18 . The idea of quantum spin Hall insulators has been extended to three-dimensional (3D) topological insulators with protected gapless surface states Fu07a ; Moore07 ; Roy09 ; Guo09 ; Weeks10 . Previous studies on 3D topological insulators have revealed that they show richer properties than those in 2D. For example, magnetoelectric properties of 3D topological insulators are described by axion electrodynamics Qi08 and, as a consequence, an image magnetic monopole is induced by an electric charge brought near the surface Qi09 .

As well as electrons, bosonic particles such as magnons Fujimoto09 ; Katsura10 ; Matsumoto11a ; Matsumoto14 ; Shindou13a ; Shindou13b ; Kim16 ; Onose10 ; Ideue12 ; Chisnell15 ; Hirschberger15 ; Han_Lee17 ; Murakami_Okamoto17 ; Kawano19 , photons Onoda04 ; Hosten08 ; Raghu08 ; Haldane08 ; Wang09a ; Ben-Abdallah16 , phonons Strohm05 ; Sheng06 ; Inyushkin07 ; Kagan08 ; Wang09b ; Zhang10 ; Qin12 ; Mori14 ; Huber16 ; Sugii17 , and triplons Rumhanyi15 ; Joshi18 , have been shown to exhibit fascinating phenomena which originate from the topology of band structures. Among them, magnons are of particular interest. The research on topological magnon systems was triggered by the prediction of the thermal Hall effect of magnons Katsura10 , which was soon observed in the pyrochlore ferromagnet Lu2V2O7 Onose10 . Other notable examples include the magnonic analog of spin Hall insulators Zyuzin16 ; Nakata17 , and topological semimetal phases as well Fransson16 ; Owerre17 ; Li16 ; Mook16 ; Su17 . The topological properties of the magnon systems with a bulk band gap are characterized by the presence of gapless surface states. In the previous study, we have defined the ${\mathbb{Z}}_{2}$ topological invariant that determines the presence or absence of magnonic helical edge states Kondo19 . The ${\mathbb{Z}}_{2}$ topological invariant is defined based on the symmetry of the pseudo-time-reversal operator whose square is $-1$. The studies of ${\mathbb{Z}}_{2}$ topological phases of magnons have so far been limited to 2D systems. However, as in the case of electrons, we expect that topologically protected surface states of magnons in 3D exhibit novel phenomena arising from the coupling of their spin and orbital motion, which never appear in 2D.

In this paper, we propose three models of magnonic analog of 3D topological insulator, which we dub the 3D topological magnon systems. They are described by a bosonic Bogoliubov-de Gennes (BdG) Hamiltonian. Such systems are intrinsically non-Hermitian, as their Hamiltonians must be diagonalized by paraunitary matrices. According to the recent classification of topological phases of non-Hermitian systems Kawabata18 ; Zhou18 , our models belong to the symmetry class AII in 3D. Following our previous work Kondo19 , we define a set of topological invariants for these systems based on Kramers pairs in bosonic systems. We discuss the correspondence of these topological invariants with the parity of the number of surface Dirac cones. One of the three models is reminiscent of the Fu-Kane-Mele model Fu07a . Since this model has inversion symmetry, one can apply a bosonic counterpart of the simplified formula for the topological invariant Fu07b , allowing us to compute it analytically. At the end of this paper, we predict that the thermal Hall effect occurs in the presence of the surface magnetic field.

The rest of the paper is organized as follows. In Sec. II, we define a set of topological invariants that characterizes the topological phases of 3D topological magnon systems and discuss the correspondence between topological invariants and the surface states. In Sec. III, we propose the first model of 3D topological magnon systems. The above correspondence is confirmed by calculating the band structure of the system with open boundary conditions. Next, in Sec. IV, we discuss a possible realization of the thermal Hall effect in the presence of a magnetic field on the surface in proximity to a normal ferromagnet. Concluding remarks are presented in Sec. V. In App. A, we show the details of the interactions in the models which are proposed in Sec. III. In App. B, the simplified expression of topological invariants for the systems with inversion symmetry is applied to this model. In App. C, we calculate analytically the energy spectrum of this model in a certain case. In App. D and E, we propose the second and third model of 3D topological magnon systems and confirm the correspondence between topological invariants and the surface states. In App. F, we show the details of Dzyaloshinskii-Moriya (DM) interaction in the third model.

II DEFINITION OF TOPOLOGICAL INVARIANT

We consider a system of non-interacting bosons in 3D lattices. Assuming translational invariance, a generic quadratic Hamiltonian describing the system is given by

[TABLE]

Here, $\bm{\beta}^{\dagger}(\bm{k})=[\beta_{1}^{\dagger}(\bm{k}),\cdots,\beta_{\mathscr{N}}^{\dagger}(\bm{k})]$ denotes boson creation operators with momentum $\bm{k}=(k_{x},k_{y},k_{z})$. The subscript $\mathscr{N}$ is a number of internal degrees of freedom in a unit cell. The $2\mathscr{N}\times 2\mathscr{N}$ Hermitian matrix $H(\bm{k})$ is a bosonic BdG Hamiltonian.

We consider the pseudo-time-reversal operator $\Theta^{\prime}$ which is introduced in our previous work Kondo19 . This operator is written as $\Theta^{\prime}=PK$, where $P$ and $K$ are a paraunitary matrix satisfying $P^{\dagger}\Sigma_{z}P=\Sigma_{z}$ and the complex conjugate operator, respectively. Here, $\Sigma_{z}$ is defined as a tensor product $\Sigma_{z}:=\sigma_{z}\otimes 1_{\mathscr{N}}$ where $\sigma_{a}$ $(a=x,y,z)$ is the $a$-component of the Pauli matrix acting on the particle-hole space and $1_{\mathscr{N}}$ is the $\mathscr{N}\times\mathscr{N}$ identity matrix. As shown in Ref Kondo19 , a $2\mathscr{N}$-dimensional complex vector $\bm{\psi}$ and $\Theta^{\prime}\bm{\psi}$ are orthogonal, i.e., $\langle\!\langle\bm{\psi},\Theta^{\prime}\bm{\psi}\rangle\!\rangle=0$ under the scalar product defined by

[TABLE]

where, $\bm{\phi}^{\dagger}$ is the adjoint of the vector $\bm{\phi}$. Then the existence of bosonic Kramers degeneracy follows from the pseudo-time-reversal symmetry:

[TABLE]

Assuming the pseudo-time-reversal symmetry, we define the ${\mathbb{Z}}_{2}$ topological invariants for 3D magnon topological systems as follows:

[TABLE]

where $n$ is a band index and $i=x$, $y$, $z$. Here, $j$ and $k$ represent two of $x$, $y$, $z$ which are different from $i$. The index $\sigma=\pm$ indicates the particle and hole spaces. The notation ${\rm EBZ}_{i,0}$ $({\rm EBZ}_{i,\pi})$ stands for the effective Brillouin zone in the $k_{i}=0$ $(k_{i}=\pi)$ plane, which is defined as $k_{i}=0$ $(k_{i}=\pi),k_{j}\in[-\pi,\pi],k_{k}\in[0,\pi]$. Its boundary is denoted by $\partial{\rm EBZ}_{i,0}$ $(\partial{\rm EBZ}_{i,\pi})$. 3D Brillouin zone and one of the effective Brllouin zones, ${\rm EBZ}_{x,0}$, are shown in Fig. 1. Here we note that any 3D lattice can be deformed so that the Brillouin zone becomes a cube as shown in Fig. 1(a) without changing the topology of the band structure. Using an index $l(=1,2)$ for Kramers doublet, $\bm{A}_{n\sigma}(\bm{k})$ and $\bm{\Omega}_{n\sigma}(\bm{k})$ are defined as

[TABLE]

where the Berry connection $\bm{A}_{nl\sigma}(\bm{k})$ and Berry curvature $\bm{\Omega}_{nl\sigma}(\bm{k})$ for bosons are given by

[TABLE]

Here, $\bm{\Psi}_{nl\sigma}(\bm{k})$ is an eigenvector of $\Sigma_{z}H(\bm{k})$. As mentioned in Ref. Kondo19 , the topological invariants of the particle and hole have the same value: $\nu_{i,0}^{n+}=\nu_{i,0}^{n-},\nu_{i,\pi}^{n+}=\nu_{i,\pi}^{n-}$. Thus, it suffices to consider the case $\sigma=+$. In the following, we drop $\sigma$ and simply write $\nu_{i,0}^{n+}$ and $\nu_{i,\pi}^{n+}$ as $\nu_{i,0}^{n}$ and $\nu_{i,\pi}^{n}$, respectively.

The number of the surface Dirac cones appearing in a bulk band gap can be understood in terms of the $\mathbb{Z}_{2}$ topological invariants. For a given bulk gap, oddness or evenness of the number of in-gap Dirac points at time-reversal-invariant momenta (TRIM) in ${\rm EBZ}_{i,0}$ $({\rm EBZ}_{i,\pi})$ is given by $\nu_{i,0}$ $(\nu_{i,\pi})$ with

[TABLE]

where $E_{n}({\bm{k}})$ is the energy of the $n$th band and $\epsilon_{0}$ is a fictitious Fermi energy that can be arbitrary as long as it lies in the gap. Since $\nu_{i,0}+\nu_{i,\pi}$ corresponds to the parity of the total number of the surface Dirac points, one has

[TABLE]

Due to this equation, only four of the six topological invariants are independent. Thus the topological phases of the system are completely characterized by the set of four topological invariants: $(\nu_{0};\nu_{x},\nu_{y},\nu_{z})$, where

[TABLE]

Following the discussion in Ref. Fu07a , for $\nu_{0}=1$, there exist odd number of surface Dirac cones in total. Such a topological phase is robust against disorder respecting pseudo-time-reversal symmetry and corresponds to the strong topological phase of electronic topological insulators. If the four topological invariants are all zero, the system is in the trivial phase. When $\nu_{0}=0$ but at least one of $\nu_{i}$ $(i=x,y,z)$ is nonzero, there exist even number of surface Dirac cones in total. However, this phase is not robust even against perturbations preserving pseudo-time-reversal symmetry, because even number of surface Dirac cones can be paired up and removed.

III MODELS

A Hamiltonian

Now we propose a model Hamiltonian of 3D topological magnon systems. We assume that the bosons from up and down spins have the same number of internal degrees of freedom. For such systems, the set of boson creation operators $\bm{\beta}^{\dagger}(\bm{k})$ in Eq. (3) takes the form:

[TABLE]

Here, the creation operators of boson from up spins $\bm{b}_{\uparrow}^{\dagger}(\bm{k})$ and those from down spins $\bm{b}_{\downarrow}^{\dagger}(\bm{k})$ are given by

[TABLE]

where $N$ is the number of the internal degrees of freedom. The operator $b_{s,i}^{\dagger}(\bm{k})$ creates a boson with spin $s$ ($s=\uparrow,\downarrow$) at site $i$. In this setup, we take the pseudo-time-reversal operator to be of the form

[TABLE]

Here we introduce our first example, which is a system in which two spins with opposite directions are localized at each site of the diamond lattice. The system is depicted in Fig. 2. The Hamiltonian of the system is written as

[TABLE]

Here, $H_{\rm DM}$ is the DM interaction between next-nearest-neighbor spins which point in the same direction. The second term $H_{J^{\prime}}$ is the antiferromagnetic interaction between two spins on the same lattice site. The third term $H_{J}$ is the bond-dependent ferromagnetic interaction between nearest-neighbor spins aligned in the same direction. The magnon band structure opens a gap due to the bond dependence of this interaction. The terms $H_{\rm XY}$ and $H_{\rm\Gamma}$ are the anisotropic XY and ${\rm\Gamma}$ interactions between next-nearest-neighbor spins pointing in opposite directions. The last term $H_{\kappa}$ is the single ion anisotropy. The interactions above realize the Fu-Kane-Mele model  Fu07a in magnon systems. Ferromagnetic Heisenberg interaction between nearest neighbor sites ($H_{J}$) plays the role of the hopping term $t\sum_{\langle i,j\rangle}c_{i}^{\dagger}c_{j}$ in Eq. (4) in Ref. Fu07a . The DM interaction ($H_{{\rm DM}}$), the Gamma interaction ($H_{{\rm\Gamma}}$), and the Heisenberg interaction with XY anisotropy ($H_{{\rm XY}}$) are considered to mimic the spin-orbit interaction $i(8\lambda_{SO}/a^{2})\sum_{\langle\!\langle i,j\rangle\!\rangle}c_{i}^{\dagger}\bm{s}\cdot(\bm{d}_{ij}^{1}\times\bm{d}_{ij}^{2})c_{j}$ in Eq. (4) in Ref. Fu07a . To sum up, the correspondence between the magnetic and the fermionic Hamiltonians is as follows:

[TABLE]

As shown in Section III C, the interaction $H_{J^{\prime}}$ does not change the topology of the band structure. The easy axis anisotropy $H_{\kappa}$ is needed to make the spin configuration in Fig. 2 stable. The details of these interactions are shown in Appendix A.

By applying the Holstein-Primakoff and Fourier transformation, we can write down the Hamiltonian in the same form as Eq. (3) with $N=2$. The matrix $H(\bm{k})$ here is proportional to the length of spins $S$. In the following, we take $S=1$ for simplicity. The explicit expression for the matrix $H(\bm{k})$ is written as follows:

[TABLE]

where $h(\bm{k})$ and $\Delta(\bm{k})$ are $4\times 4$ matrices defined as

[TABLE]

with

[TABLE]

These $\Gamma$ matrices satisfy the standard anticommutation relations: $\{\Gamma^{(i)},\Gamma^{(j)}\}=2\delta_{ij}1_{4}$. The six components in Eq. (27) are defined as

[TABLE]

where $k_{i}=\bm{k}\cdot\bm{a}_{i}$ and $k_{ij}=k_{i}-k_{j}$ and all the parameters are real. The vectors $\bm{a}_{i}$ $(i=1,2,3)$ are the three lattice prime vectors of the diamond lattice. We note that the Hamiltonian (26) satisfies the pseudo-time-reversal symmetry (5).

B $\mathbb{Z}_{2}$ topological invariants

Let us now compute the $\mathbb{Z}_{2}$ topological invariants. Thanks to the inversion symmetry of the system, one can compute them analytically by using a simplified formula which can be thought of as the bosonic counterpart of the formula derived in Ref. Fu07b (see Appendix B for details). The results are summarized in Table 1. We have checked that the topological invariants obtained are the same as those obtained by evaluating Eq. (6) numerically. Table 1 suggests that the system is in the strong topological phase, i.e., odd number of Dirac cones exist between the top and bottom bands in the particle space (the correspondence of topological invariants to the number of Dirac cones is confirmed later). The bulk band structure with the same parameters as those of Table 1 is shown in Fig. 3(a). Since the system has both the pseudo-time-reversal symmetry and inversion symmetry, the energy spectrum is doubly degenerate over the whole Brillouin zone. We note that nontrivial topological phases can be realized when at least two of $J_{-},D$, and $\Gamma$ are nonzero Energy_analytical .

C Phase diagram and magnon band structure

We now construct a phase diagram of the diamond lattice system by computing the simplified formula of topological invariants analytically. The resulting phase diagram is shown in Fig. 3(b) (see Appendix B for a detailed derivation). By changing the coupling constants, three phases: strong topological, weak topological, and trivial phases are all realized. The band structures for a slab with (100) face of the system which is deformed to a cubic lattice are shown in Fig. 4. As we expect from the general arguments of the relation between the topological invariants and the number of the surface Dirac points. In weak topological phases $(0;111)$ and $(0;100)$ (Fig. 4(a) and (b), respectively), there are even number (2 and 0, respectively) of surface Dirac cones. On the other hand, in strong topological phases $(1;111)$ and $(1;100)$ (Fig. 4(c) and (d), respectively), there are odd number (1 and 3, respectively) of surface Dirac cones. On general grounds, it is expected that the surface states arising in the strong topological phase exhibit magnon spin-momentum locking Okuma17 even in the collinear magnet due to the terms $H_{\rm\Gamma}$ and $H_{\rm XY}$ which break spin conservation.

Let us now show that the magnonic system with nonzero $J^{\prime}$ is in the topological phase corresponding to that of the Fu-Kane-Mele model Fu07a , whose Hamiltonian is essentially the same as Eq. (26) with $J^{\prime}=0$. Numerical results suggest that the band gap of $H({\bm{k}})$ closes only at TRIM, in which case this magnon system is adiabatically connected to the system with $J^{\prime}=0$ without closing a gap. This can be seen as follows. We first note that the band gap at TRIM $\bm{k}=\bm{\Gamma}$ can be written as

[TABLE]

(see Appendix B for a detailed derivation). Clearly, $\Delta E=0$ if $d_{1}(\bm{\Gamma})=0$ and $\Delta E>0$ otherwise. Now suppose that the gap opens at $J^{\prime}\neq 0$. This implies $d_{1}(\bm{\Gamma})\neq 0$ at all TRIM. Then one finds that the gap does not close even if one takes $J^{\prime}\to 0$. This proves the desired result.

D Other models

So far we have focused on the model on the diamond lattice. However, the physics of the 3D topological magnon systems is not limited to this lattice. To illustrate this, we have constructed two other examples: one is defined on the perovskite lattice and the other on the the pyrochlore lattice. In both cases, the validity of the relation between the set of topological invariants and the number of the surface Dirac points has been confirmed (see Appendix D and E for details).

IV. THERMAL HALL EFFECT

We now turn to discuss the physical implications of surface Dirac cones in the strong topological phase. For electronic systems, previous studies have shown that the surface Dirac dispersions in a 3D topological insulator can be gapped out by applying a magnetic field on the surface Sinitsyn07 ; Lu10 . In general, such states have nonzero Berry curvature, leading to the surface quantum Hall effect. Here we show that the magnonic counterpart of this phenomenon occurs in our system. To this end, we consider a heterostructure of a ferromagnet and a 3D topological magnon system shown in Fig. 5. To be concrete, we take the diamond lattice system as a 3D topological magnon system. We assume that the magnetic field emanating from the ferromagnet is weak enough not to affect the magnetic order in the substrate.

First, we derive the effective Hamiltonian of the surface of the system. Let us consider the projection onto the subspace of two (particle) bands above and below the Dirac point we focus on. Denoting the indices of the two bands as $n$ and $n+1$, we define the set of the wave functions of these bands as

[TABLE]

The Hamiltonian projected onto this subspace is a $4\times 4$ matrix and written as

[TABLE]

Expanding the Hamiltonian $H_{n}(\bm{k})$ to the second order in $k_{z}$ and the first order in $(k_{x}-\pi)$ and ($k_{y}-\pi)$, we denote the terms depending on $(k_{x}-\pi)$ or $(k_{y}-\pi)$ as $H_{xy}(k_{x},k_{y})$. The other terms are denoted by $H_{z}(k_{z})$. We write the energy of the center of the Dirac cone as $E_{0}$. The state $\bm{\psi}(z)$ which is localized on the $(001)$ surface is obtained as the eigenvector of the Hamiltonian $H_{z}(-{\rm i}\partial_{z})$ with the eigenvalue $E_{0}$:

[TABLE]

Here we take $\bm{\psi}(z)$ as

[TABLE]

The parameter $\lambda$ is obtained by solving

[TABLE]

This is an eighth order equation in $\lambda$ and has eight solutions. Four of them are positive and the other are negative. The negative solutions correspond to states localized on the $(00\bar{1})$ surface. Here we focus on the positive solutions corresponding to states localized on the $(001)$ surface. Due to the time-reversal symmetry, each two of the eight solutions are the same. Thus, we write two positive solutions as $\lambda_{1}$ and $\lambda_{2}$ and define the corresponding two constant vectors as $\bm{\psi}_{1}$ and $\bm{\psi}_{2}$. The pseudo-time-reversal operator in the subspace of $n$th and $(n+1)$th bands is defined as $\Theta=({\rm i}\sigma_{y}\otimes 1_{2})K$. The vectors $\Theta\bm{\psi}_{1}$ and $\Theta\bm{\psi}_{2}$ are also the solutions with $\lambda_{1}$ and $\lambda_{2}$, respectively. Then the wave function which satisfies Eq. (31) is generally written as

[TABLE]

Since $\lambda_{1,2}>0$, this satisfies the condition $\bm{\psi}(-\infty)=0$. Using the other boundary condition that the wave function vanishes on the surface $z=0$, i.e. $\bm{\psi}(0)=0$, we obtain four simultaneous equations for the coefficient $\alpha_{1},\beta_{1},\alpha_{2}$ and $\beta_{2}$:

[TABLE]

Since the determinant of the matrix $\left(\bm{\psi}_{1},\Theta\bm{\psi}_{1},\bm{\psi}_{2},\Theta\bm{\psi}_{2}\right)$ is zero, the above equation has a nontrivial solution. The overall factor of the four coefficients $\alpha_{1},\beta_{1},\alpha_{2}$ and $\beta_{2}$ is determined by the normalization condition:

[TABLE]

From the coefficients obtained, we have the effective surface Hamiltonian

[TABLE]

where the matrix element $\bra{\bm{\psi}}H_{xy}(k_{x},k_{y})\ket{\bm{\psi}}$ is defined as

[TABLE]

The other three matrix elements are defined similarly. By applying the magnetic field $\bm{B}=B\bm{e}_{z}$, the additional term $-BS\sigma_{z}$ appears in the Hamiltonian. The band structure of the surface is obtained by diagonalizing the Hamiltonian $H_{\rm eff}(k_{x},k_{y})$. The Berry curvature of the system is defined as $\Omega_{n}^{z}(k_{x},k_{y})=2{\rm Im}\left[(\partial_{k_{x}}\bm{\psi}_{n}^{\dagger}(k_{x},k_{y}))(\partial_{k_{y}}\bm{\psi}_{n}(k_{x},k_{y}))\right]$. Figure 6 shows the band structure of the surface Hamiltonian without and with the surface magnetic field. We can see that the surface state can be gapped out by applying the surface magnetic field. Figure 7 shows the Berry curvature of the top and bottom bands under magnetic field, respectively. Since the surface state has nonzero Berry curvature, the thermal Hall coefficient calculated in Ref Matsumoto11b :

[TABLE]

is expected to be nonzero. Here, $k_{B},\hbar,V,$ and $T$ are the Boltzmann constant, the Planck constant, the volume of the system, and the average temperature of the system, respectively. The function $c_{2}(\rho)$ is defined as $c_{2}(\rho)=(1+\rho)\left(\log\frac{1+\rho}{\rho}\right)^{2}-\left(\log\rho\right)^{2}-2{\rm Li}_{2}(-\rho)$, where ${\rm Li}_{2}(-\rho)$ and $\rho(E_{n}(k_{x},k_{y}))$ are the dilogarithm function and the Bose distribution function, respectively.

V. SUMMARY

In this paper, we proposed three models of 3D topological magnon systems that host topologically protected surface Dirac cones. The models have pseudo-time-reversal symmetry which allows us to define the set of $\mathbb{Z}_{2}$ topological invariants that can distinguish different phases of these systems. For a model on the diamond lattice, we computed the topological invariants analytically by exploiting the inversion symmetry of the system. The analysis of these topological invariants revealed that the model exhibits strong topological, weak topological, and trivial phases, in which the number of surface Dirac points is odd, even, and zero, respectively.

We also studied the topological invariants in the other two models numerically and found the same correspondence. As a physical consequence of the presence of surface Dirac dispersions, we discussed the thermal Hall effect in the composite system consisting of a ferromagnet and a 3D topological magnon system. In this setup, the Zeeman field at the interface gaps out the surface Dirac dispersions, leading to the nonvanishing Berry curvature, and hence a nonzero thermal Hall response. Finally, we remark that the present formulation of the $\mathbb{Z}_{2}$ topological invariants can be applied to other bosonic systems such as phonons and photons, as long as they respect pseudo-time-reversal symmetry. In this respect, it would be interesting to study bosonic excitations in spin liquids or paramagnets, which are entirely different from the magnon excitations. We expect that bosonic $\mathbb{Z}_{2}$ topological phases in such systems can be identified by combining our approach with the Schwinger-boson mean-field theory Sachdev92 ; Wang06 ; Liu19 .

Acknowledgements

This work was supported by JSPS KAKENHI Grants No. JP18H04478, No. JP18K03445, JP17K14352, and JP18H04220. H. Kondo was supported by the JSPS through Program for Leading Graduate Schools (ALPS).

Appendix A: Details of the interactions

In this part, we describe the details of the interactions of the first model (19). They are written as follows:

[TABLE]

where $\bm{S}_{i,{\rm u}}$ and $\bm{S}_{i,{\rm d}}$ are the operators of up and down spins localized at the site $i$, respectively. The operator $\bm{S}_{s}(\bm{R},X)$ $(s={\rm u,d})$ is the another expression of $\bm{S}_{i,s}$ $(s={\rm u,d})$ for the site $i$ on $X$ sublattice in the unit cell labeled by the lattice vector $\bm{R}$. The vector $\bm{a}_{ij}$ is defined as $\bm{a}_{ij}=\bm{a}_{i}-\bm{a}_{j}$. The DM vectors $\bm{D}_{i}$ ($\bm{D}_{ij}$) is written as $\bm{D}_{i}=D(\bm{d}^{1}_{i}(\bm{R})\times\bm{d}^{2}_{i}(\bm{R}))/|\bm{d}^{1}_{i}(\bm{R})\times\bm{d}^{2}_{i}(\bm{R})|$ $\left(\bm{D}_{ij}=D(\bm{d}^{1}_{ij}(\bm{R})\times\bm{d}^{2}_{ij}(\bm{R}))/|\bm{d}^{1}_{ij}(\bm{R})\times\bm{d}^{2}_{ij}(\bm{R})|\right)$, where $\bm{d}^{1,2}_{i}(\bm{R})$ $(\bm{d}^{1,2}_{ij}(\bm{R}))$ are the two nearest neighbor bond vectors traversed between sites $(\bm{R},A)$ and $(\bm{R}+\bm{a}_{i},A)$ ($(\bm{R},A)$ and $(\bm{R}+\bm{a}_{ij},A)$). Figure 8 shows the details of $\bm{d}^{1,2}_{i}(\bm{R})$. The bond dependence of the Heisenberg interaction (48) originates from pressure-induced lattice distortion. Here $\bar{\bm{D}}_{i}$ and $\bar{\bm{D}}_{ij}$ are written as $\bar{\bm{D}}_{i}=\bm{D}_{i}/D$ and $\bar{\bm{D}}_{ij}=\bm{D}_{ij}/D$, respectively.

Appendix B: Simplified expressions for the topological invariants of the diamond lattice system

If a system has inversion symmetry, the expression of the topological invariants for topological insulators can be greatly simplified Fu07b . The same discussion can be applied to magnon systems. In this part, by using inversion symmetry, we analytically calculate the topological invariants for diamond lattice system. The Hamiltonian (26) satisfies the following inversion symmetry:

[TABLE]

where $R$ is an inversion operator defined as $R:=1_{2}\otimes 1_{2}\otimes\sigma_{x}$. Following the discussion in Ref. Fu07b , topological invariants for 3D topological magnon systems with inversion symmetry can be written as

[TABLE]

where $i=x,y$, and $z$. Since $\Sigma_{z}H(\bm{k})$ commutes with the inversion operator $R$ at TRIM: $\bm{\Gamma}_{m}=\pi(n_{1},n_{2},n_{3})$, an energy eigenvector $\bm{\Psi}_{n,1,+}(\bm{\Gamma}_{m})$ can be chosen to be an eigenvector of $R$. Here, we write the eigenvalue of $R$ as $\xi_{n}(\bm{\Gamma}_{m})$. Then, $\delta_{m=(n_{1}n_{2}n_{3})}$ is defined as the product of $\xi_{n}(\bm{\Gamma}_{m})$ over the bands below the fictitious Fermi energy:

[TABLE]

Now, we assume that the fictitious Fermi energy lies between two bands in the particle space. In this case, $\delta_{m=(n_{1}n_{2}n_{3})}$ is an eigenvalue of $R$ with eigenvector of the lower band $\bm{\Psi}_{2,1,+}(\bm{\Gamma}_{m})$. In the following, we obtain the expression of $\bm{\Psi}_{2,1,+}(\bm{\Gamma}_{m})$. At the TRIM, the Hamiltonian of the diamond lattice system (163) is written as

[TABLE]

The eigenvectors of the matrix $h_{1}(\bm{\Gamma})$ are written as:

[TABLE]

where the corresponding eigenvalues are $\epsilon_{\pm}(\bm{\Gamma})=d_{0}\pm|d_{1}(\bm{\Gamma})|$. Here, $\gamma(\bm{k})$ is defined as $\gamma(\bm{k})=-d_{1}(\bm{k})-{\rm i}d_{2}(\bm{k})$. The concrete expressions of $d_{0}$ and $d_{i}(\bm{k})$ $(i=1,\cdots,5)$ are given in the main text. We note that imaginary part of $\gamma(\bm{k})$, i.e. $d_{2}(\bm{k})$, vanishes at TRIM. If we write the eigenvector of $\Sigma_{z}H(\bm{\Gamma})$ as

[TABLE]

the eigenvalue equation can be written as follows:

[TABLE]

This results in the following relation:

[TABLE]

In order for the equation to have a nontrivial solution, the determinant of the matrix in the left hand side must be zero. Then, we obtain the following equation:

[TABLE]

Therefore, we obtain the energy eigenvalues in the particle space:

[TABLE]

We can see that $E_{+}(\bm{\Gamma})$ and $E_{-}(\bm{\Gamma})$ correspond to the energies of the top and bottom bands, respectively. Now we focus on the bottom band. By using $\theta(\bm{\Gamma})$ defined as

[TABLE]

coefficients $\alpha$ and $\beta$ of the normalized wave function of the bottom band are written as $\alpha=\cosh{(\theta(\bm{\Gamma}))}$ and $\beta=\sinh{(\theta(\bm{\Gamma}))}$, respectively. Then, the eigenvector of the bottom band $\bm{\Psi}_{2,1,+}(\bm{\Gamma})$ are given by:

[TABLE]

We can see that this is an eigenvector of $R$ with eigenvalue ${\rm sgn}[\gamma(\bm{\Gamma})]$. Then, $\delta_{m=(n_{1}n_{2}n_{3})}$ in Eq. (54) can be written as

[TABLE]

By using this, the strong index for the diamond lattice system is given by

[TABLE]

The other three indices are written as follows:

[TABLE]

Appendix C: The analytical expressions of energy eigenvalues

In the case of $\Gamma=0$, the energy spectrum of the Hamiltonian of the diamond lattice system (19) can be obtained analytically. By using the unitary matrix $U(\bm{k})$ (details of $U(\bm{k})$ are shown later), $\Sigma_{z}H(\bm{k})$ can be written as follows:

[TABLE]

where $Q_{n}(\bm{k})$ $(n=1,2)$ are defined by

[TABLE]

Here $\lambda_{n,\pm}(\bm{k})$ is written as

[TABLE]

We note that the matrix $U^{\dagger}(\bm{k})(\Sigma_{z}H(\bm{k}))^{2}U(\bm{k})$, whose eigenvalues are equal to the square of the ones of $\Sigma_{z}H(\bm{k})$, is the block diagonal:

[TABLE]

This suggests that the square of the energy eigenvalue $E^{2}(\bm{k})$ can be obtained as the eigenvalues of $Q_{1}(\bm{k})Q_{2}(\bm{k})$:

[TABLE]

Since this matrix can be divided into two parts, the eigenvalues are calculated from the following equation:

[TABLE]

Then, one has

[TABLE]

Finally, the eigenvalues of $\Sigma_{z}H(\bm{k})$ are written as follows:

[TABLE]

where $\rho_{1},\rho_{2},\rho_{3}=\pm$. We note that due to the pseudo-time-reversal and inversion symmetry, these eigenvalues are doubly degenerate and satisfy $E_{\rho_{1},+,\rho_{3}}(\bm{k})=E_{\rho_{1},-,\rho_{3}}(\bm{k})$.

In the following, we show how to construct the matrix $U(\bm{k})$. This is written as the product of three unitary matrices $U_{1},U_{2},$ and $U_{3}(\bm{k})$, i.e. $U(\bm{k})=U_{1}U_{2}U_{3}(\bm{k})$. The first matrix $U_{1}$ is the one which diagonalizes $C:=\sigma_{y}\otimes\sigma_{y}\otimes 1_{2}$. Since $C$ anticommutes with $\Sigma_{z}H(\bm{k})$, $U_{1}$ makes $\Sigma_{z}H(\bm{k})$ off diagonal:

[TABLE]

where

[TABLE]

[TABLE]

The second matrix is given by $U_{2}=1_{2}\otimes u$. Here $u$ is a $4\times 4$ matrix and diagonalizes $c:=\sigma_{x}\otimes\sigma_{z}$ which anticommutes with $d_{0}1_{4}+R_{1,2}(\bm{k})$. By the unitary transformation using $u$, $d_{0}1_{4}+R_{1,2}(\bm{k})$ becomes off diagonal:

[TABLE]

[TABLE]

where

[TABLE]

[TABLE]

Since $r_{ij}(\bm{k})$ $(i,j=1,2)$ commute with each other, they are diagonalized by the same unitary matrix $v(\bm{k})$ as follows:

[TABLE]

By using $v(\bm{k})$, the third matrix $U_{3}(\bm{k})$ is defined as $U_{3}(\bm{k})=1_{4}\otimes v(\bm{k})$. Now we can construct the matrix $U(\bm{k})$ and obtain the unitary-transformed Hamiltonian Eq. (103).

Appendix D: Second model of 3D topological magnon systems: perovskite lattice system

The second example of 3D magnon topological systems is a system in which two spins with opposite directions are localized at each site of perovskite lattice. The system is depicted in Fig. 9. The Hamiltonian of the system are written as

[TABLE]

Here, $H_{\rm DM},H_{J}$ are the DM and ferromagnetic Heisenberg interactions between nearest-neighbor spins aligned in the same direction. The third term $H_{\rm XY}$ is the anisotropic XY interaction between nearest-neighbor spins pointing in opposite directions. The other terms $H_{J^{\prime}}$ and $H_{\kappa}$ are the antiferromagnetic interaction between spins at the same lattice site and the single ion anisotropy, respectively. We note that the single-ion anisotropy is introduced to open a band gap.

The explicit forms of the terms are given by

[TABLE]

Here, we take the directions of ordered spins as $\pm\bm{s}=\pm(s_{1},s_{2},s_{3})$ and write corresponding spin operators at site $i$ as $\bm{S}_{i,{\rm u}},\bm{S}_{i,{\rm d}}$, respectively. The Hamiltonian of anisotropic XY interaction is written as follows:

[TABLE]

The diagonal matrices $J_{\rm XY}$ and $J_{\rm YX}$ are $3\times 3$ matrices and written as $J_{\rm XY}={\rm diag}[J_{x},J_{y},0]$ and $J_{\rm YX}={\rm diag}[J_{y},J_{x},0]$, respectively. Assuming the ligands are located on the center of the octahedra of the perovskite lattice, we difine the DM vector $\bm{D}_{ij}$ for the nearest neighbor sites $i,j$ as $\bm{D}_{ij}=D(\bm{d}_{i}\times\bm{d}_{j})/|\bm{d}_{i}\times\bm{d}_{j}|$, where $\bm{d}_{i}$ and $\bm{d}_{j}$ are the vectors from sites $i,j$ to the nearest ligand.

Here, the most general Hamiltonian satisfying Eq. (5) takes the form:

[TABLE]

where $h_{i}(\bm{k})$ and $\Delta_{i}(\bm{k})$ for $i=1,2$ are $3\times 3$ matrices and satisfy $h_{1}^{\dagger}(\bm{k})=h_{1}(\bm{k}),\Delta_{1}^{\dagger}(\bm{k})=\Delta_{1}(\bm{k}),h_{2}^{T}(\bm{k})=-h_{2}(-\bm{k}),$ and $\Delta_{2}^{T}(\bm{k})=\Delta_{2}(-\bm{k})$. Applying the Holstein-Primakoff transformation, the Hamiltonian (154) boils down to the same form as Eq. (163). We take the operators $\bm{b}_{\uparrow}^{\dagger}(\bm{k})$ and $\bm{b}_{\downarrow}^{\dagger}(\bm{k})$ in Eq. (17) as

[TABLE]

Writing the length of the spin as $S$, the concrete expressions of $h_{1}(\bm{k}),h_{2}(\bm{k}),\Delta_{1}(\bm{k}),$ and $\Delta_{2}(\bm{k})$ are given by

[TABLE]

where

[TABLE]

Here $J_{\pm}$ is defined as $J_{\pm}=(J_{x}\pm J_{y})/2$.

Using the numerical implementation described in Ref. Fukui07 , we calculate the topological invariants of the system. The results are shown in Table 2. We note that the realization of nontrivial topological phases depend on the direction of ordered spins $\bm{s}$. For example, in the case of $\bm{s}=(0,0,1)$, the band structure is topologically trivial, while the index $\nu_{0}^{n}=1$ $(n=1,2)$ in Table 2 implies that the system is in the strong topological phase if $s_{1}=-(2+\sqrt{2})/4,s_{2}=(2-\sqrt{2})/4,s_{3}=1/2$. In order to see the correspondence between the topological invariants in Table 2 and the number of Dirac cones, let us plot the band structure of the system. The band structure of the bulk system and a slab with (001) face are shown in Fig. 10. As far as we investigate this system, the strong topological phase is realized only in the system with the indirect gap.

As shown in Fig. 10(b), a single Dirac cone exists at the ${\rm M}_{3}$ point of the upper gap. This corresponds to the summation of the strong indices of the lower and middle bands $\nu_{0}^{2}+\nu_{0}^{3}=1$. This Dirac cone is robust against disorder which does not break the pseudo-time-reversal symmetry. The topological invariants of the lower band are all zero. This indicates the absence of the Dirac cones in the gap between the lower and middle bands.

Appendix E: Third model of 3D topological magnon systems: pyrochlore lattice system

We give another model of 3D magnon topological systems. This is the system in which the same number of up and down spins are localized at the same site of the pyrochlore lattice. The schematic picture of the system is shown in Fig. 11.

The Hamiltonian of the system is written as

[TABLE]

where $H_{\rm DM}$ and $H_{J}$ are DM and ferromagnetic interactions between nearest neighbor spins which point in the same direction. The third term $H_{\rm XY}$ is the anisotropic XY interaction between nearest neighbor spins pointing in opposite directions. The details of the DM interactions in the pyrochlore lattice are shown in Appendix F. The other terms $H_{J^{\prime}}$ and $H_{\kappa}$ are the antiferromagnetic interaction between two spins on the same site and the single ion anisotropy, respectively. The single ion anisotropy opens the gap of the band structure. The Hamiltonians $H_{\rm DM},H_{J},H_{J^{\prime}}$, and $H_{\kappa}$ are given by

[TABLE]

where operators $\bm{S}_{i,{\rm u}}$ and $\bm{S}_{i,{\rm d}}$ denote up and down spins on site $i$, respectively. The anisotropic XY interaction $H_{\rm XY}$ is written as

[TABLE]

Applying the Holstein-Primakoff transformation, we rewrite the Hamiltonian (185) in the same form as Eq. (163). We take the operators $\bm{b}_{\uparrow}^{\dagger}(\bm{k})$ and $\bm{b}_{\downarrow}^{\dagger}(\bm{k})$ in Eq. (17) as

[TABLE]

The concrete expressions of $h_{1}(\bm{k}),h_{2}(\bm{k}),\Delta_{1}(\bm{k})$, and $\Delta_{2}(\bm{k})$ are given by

[TABLE]

where $\gamma_{i}(\bm{k})=1+e^{{\rm i}k_{i}},\gamma_{ij}(\bm{k})=1+e^{{\rm i}k_{ij}}$, and $\tilde{D}=\frac{\sqrt{3}}{6}D$.

Using the numerical implementation described in Ref. Fukui07 , we calculate the topological invariants of the system. The results are shown in Table 3. Table 3 suggests that the system has even, even, and odd number of Dirac cones in the highest, second highest, and the lowest band gaps, respectively. In order to confirm this, we calculate the band structure of the system. The band structure of the bulk system and a slab with (001) face is shown in Fig. 12. A single Dirac cone exists at the ${\rm M_{2}}$ point in the lowest band gap. This corresponds to the strong topological index $\nu_{0}^{4}=1$. The second highest and the lowest band gaps have two Dirac cones. This means that the summations of the strong topological index over the lowest two and three bands are zero:

[TABLE]

Appendix F: Details of the DM interaction in the pyrochlore lattice system

Here we write the details of DM vectors in the pyrochlore lattice system. The nearest neighbor bonds of the pyrochlore lattice are shown in Fig. 13. Here $X(\bm{R})$ is the $X$ sublattice in the unit cell with the lattice vector $\bm{R}$.

We write the DM vector between two sites $X(\bm{R})$ and $X^{\prime}(\bm{R^{\prime}})$ as $\bm{D}(X(\bm{R}),X^{\prime}(\bm{R^{\prime}}))$($\equiv\bm{D}_{ij}$). The unit vector from the site $X(\bm{R})$ to $X^{\prime}(\bm{R})$ is defined as $\overrightarrow{XX^{\prime}}$. Each DM vector of the pyrochlore lattice is written as

[TABLE]