Orbital Edelstein effect in topological insulators
OOrbital Edelstein effect in topological insulators
Ken Osumi, Tiantian Zhang,
1, 2 and Shuichi Murakami
1, 2, ∗ Department of Physics, Tokyo Institute of Technology, Tokyo 152-8551, Japan TIES, Tokyo Institute of Technology, Tokyo 152-8551, Japan (Dated: February 4, 2021)We theoretically propose a gigantic orbital Edelstein effect in topological insulators and interpretthe results in terms of topological surface currents. We numerically calculate the orbital Edelsteineffect for a model of a three-dimensional Chern insulator as an example. Furthermore, we calcu-late the orbital Edelstein effect as a surface quantity using a surface Hamiltonian of a topologicalinsulator, and numerically show that it well describes the results by direct numerical calculation.We find that the orbital Edelstein effect depends on the local crystal structure of the surface, whichshows that the orbital Edelstein effect cannot be defined as a bulk quantity. We propose that Cherninsulators and Z topological insulators can be a platform with a large orbital Edelstein effect be-cause current flows only along the surface. We also propose candidate topological insulators forthis effect. As a result, the orbital magnetization as a response to the current is much larger intopological insulators than that in metals by many orders of magnitude. In recent years, new responses leading to orbital mag-netization have been proposed in systems without in-version symmetry . One of the focuses is conver-sion of electron current and magnetization on crys-tal structure with low symmetry. Among such pro-posals are orbital Edelstein effect (OEE) ,i.e. current-induced orbital magnetization, and the gy-rotropic effect . These effects have similar response co-efficients. In particular, OEE is an orbital analog of theEdelstein effect . OEE emerges even in systems with-out spin-orbit interactions . In particular, the OEEemerges in crystals with a chiral structure , similar tothe phenomenon in which the solenoid creates a mag-netic field when a current flows. As a similar context,the recent finding of spin-selective electron transportthrough chiral molecules, the so-called chirality-inducedspin selectivity (CISS) effect , suggests an alterna-tive method of using organic materials as spin filters forspintronics applications.In the OEE, the electric field induces the magneti-zation, which may look similar to the magnetoelectriceffect . Nonetheless, in the spin and OEE, metal-lic systems are considered, and nonequilibrium electrondistribution by the electric field is a key to generatemagnetization. In this sense, the Edelstein effect canbe called kinetic magnetoelectric effect. On the otherhand, in the magnetoelectric effect, insulating systemsare considered and the system stays in equilibrium evenunder the electric field or the magnetic field. Corre-spondingly, symmetry requirements are different. In themagnetoelectric effect, inversion and time-reversal sym-metries should be broken but their product symmetry ispreserved. On the other hand, in the Edelstein effect, theinversion symmetry should be broken. It is seen in chiralsystems , and in polar systems .Magnetoelectric tensors may require careful considera-tion of boundaries. While orbital magnetization is inde-pendent on the boundary , the orbital magnetizationwhen an electric field is applied may not have such prop-erties. The general orbital magnetoelectric response depends on the boundary . Therefore, it is important to study the effect of the boundary of the response oforbital magnetization.In this paper, we investigate OEE in topological in-sulators such as three-dimensional Chern insulators and Z topological insulators ( Z -TIs) in which the cur-rents are localized on the surface. First, we calculate theOEE in three-dimensional topological insulators with chi-ral crystal structure. Second, we derive the OEE basedon the surface Hamiltonian, and we show that this effectdepends on surface states. Finally, we propose candi-date materials for this effect and estimate the values ofthe OEE. By comparing the results with the results in achiral semiconductor tellurium , we show that in topo-logical insulators the orbital magnetization as a responseto the current is much larger than metals by many ordersof magnitude. ResultsFormulation for OEE.
We consider a crystal in a shapeof a cylinder along the z -axis, and calculate its orbitalmagnetization along the z -axis generated by the currentalong the z -axis. Let c be the lattice constant along the z -axis. We introduce the velocity operator v as v = − i ¯ h [ r , H ] , where r is the position operator and H is theHamiltonian. In the limit of the system length along the z -axis to be infinity, the orbital magnetization at zerotemperature is M z = 12 π Z π/c − π/c dk z S N X n f ( E n ( k z )) × (cid:16) − e (cid:17) h ψ n ( k z ) | ( r × v ) z | ψ n ( k z ) i , (1)where − e is the electron charge, | ψ n ( k z ) i and E n arethe n th occupied eigenstates and energy eigenvalues of H at the Bloch wavenumber k z , respectively. f ( E ) isthe distribution function at the energy E , S is the crosssection of the crystal along the xy -plane and N is thenumber of occupied states.Then, within the Boltzman approximation, the appliedelectric field E z changes f ( E ) from f ( E ) into f ( E ) = f ( E )+ eτE z ¯ h ∂f ( E ) ∂k z in a linear order in E z , where τ is the a r X i v : . [ c ond - m a t . m e s - h a ll ] F e b a edc b z y x k z k y k x RXMΓZ A y x L x L y L x =10a, L y =30aL x =20a, L y =30aL x =30a, L y =30a μ/t y μ/t y L x =10a, L y =30aL x =20a, L y =30aL x =30a, L y =30a hgf FIG. 1:
The model of a Chern insulator with a chiral structure and model calculation of the OEE. ( a )Individual layer of the model forming a square lattice. The blue regions surrounded by the broken line are the unitcells consist of four sublattices. ( b ) Schematic picture of the chiral hopping (red) between the two neighboringlayers. These hoppings form structure similar to right-handed solenoids. ( c ) Brillouin zone of our model withhigh-symmetry points. ( d, e ) Energy bands for the Hamiltonian H with parameters ( d ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y and ( e ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y . ( f ) One-dimensional modelwith periodic boundary condition in z direction. In order to see the boundary effect on OEE, the outermost layerson the xz surface has no chiral hoppings and those on the yz surface has chiral hoppings. ( g, h ) OEE calculatedwith parameters ( g ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y and ( h ) t x = b x = t y , t y = b y = m , t = 0 . t y , and t = 0 . t y .relaxation time assumed to be constant and f ( E ) is theFermi distribution function f ( E ) = ( e β ( E − µ ) + 1) − ,β = 1 /k B T, k B is the Boltzman constant and µ is thechemical potential. Then the orbital magnetization isgenerated as M OEE z = 12 π Z π/c − π/c dk z S N X n eτ E z ¯ h ∂f ( E n ( k z )) ∂k z × (cid:16) − e (cid:17) h ψ n ( k z ) | ( r × v ) z | ψ n ( k z ) i . (2)This is the OEE. This calculation method is differentfrom that in the previous study on metals , where bulkcontribution in a system infinite along x and y directionsare calculated. Model calculation on a Chern insulator.
As an ex-ample of a topological insulator, we consider an orbitalmagnetization in a Chern insulator with a chiral crys-tal structure. For this purpose, we introduce a three-dimensional tight-binding model of a layered Chern in- sulator, as shown in Fig. 1a, connected via right-handedinterlayer chiral hoppings (Fig. 1b). Each layer forms asquare lattice within the xy -plane, with a lattice constant a , and they are stacked along the z -axis with a spacing c .as shown in detail in Methods. The Brillouin zone andthe band structure is shown in Figs. 1c-e. We set theFermi energy in the energy gap.We calculate OEE in a one-dimensional quadrangularprism with xz and yz surfaces shown in Fig. 1f (see Meth-ods), with its results in Figs. 1g and 1h with the interlayerhopping t = 0 . t y and t = 0 . t y , respectively, for sev-eral values of the system size, L x and L y , representingthe lengths of the crystal in the x and y directions. Thus,the OEE is affected by boundaries and system size, andthis size dependence remains even when the system sizeis much larger than the penetration depth of topologicalsurface states. Therefore, OEE cannot be defined as abulk quantity. Later, we give an interpretation on thischaracteristic size dependence. Surface theory of OEE for a slab.
In topological ...... ............y x...... ............y x
Model IModel II μ/t y Model IModel II μ/t y y=y - y=y + π/c-π/c-π/a π/a k x k z y=y - y=y + π/c-π/c-π/a π/a k x k z Model IIModel I
Model IModel II k z k z k x E valence bandconduction band E - =E - (k x k z ) y + y y + yy y y - xz v vE μ k x k z Fermi surface G x = C a b ced gf ih FIG. 2:
A Chern insulator with a chiral structure and OEE in slab systems. ( a ) Surface velocity of thetopological surface state under an electric field. ( b, c ) Band structure of the topological chiral surface states on the y = y − surface. Their ( b ) dispersion and ( c ) Fermi surface are shown. ( d ) Slab model I with chiral hopping on thesurface. ( e ) Slab model II with no chiral hopping on the surface. ( f, g ) The fermi surfaces for the slab models ( f ) Iand ( g ) II with parameters t x = t y = m = b x = b y , t = 0 . t y , t = 0 . t y and µ = 0. ( h, i ) OEE for the slab modelsI and II with parameters ( h ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y and ( i ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y .insulators such as Chern insulators, only the topologicalsurface states can carry a current. Here we calculate theOEE using an effective Hamiltonian for the crystal sur-face. Thereby, we can capture natures of OEE throughthis surface theory. We consider slab systems, with itssurfaces on y = y ± ( y + > y − ). The slab is sufficientlylong along the x and z directions and we impose peri-odic boundary conditions in these directions. To inducethe orbital magnetization M OEE z, slab , we apply an electricfield E z in the z direction. Due to the interlayer chi-ral hoppings, the surface current acquires a nonzero z -component (Fig. 2a).Let E − = E − ( k x , k z )(= E − ( k x , k z + πc )) be the surface-state dispersion on the y = y − surface as shown in Fig. 2band c. For simplicity, we assume C z symmetry of thesystem. Then the surface state dispersion on the y = y + surface is given by E + = E + ( k x , k z ) = E − ( − k x , k z ).Here, we assume that the surface states are sharply local-ized at y = y ± , namely, we ignore finite-size effects dueto a finite penetration depth. Then, we rewrite equation (2) to M OEE z, slab = e τ E z h Z π/c − π/c dk z (2 π ) × ∂E ( k x , k z ) ∂k z sgn (cid:18) ∂E ( k x , k z ) ∂k x (cid:19) (cid:12)(cid:12)(cid:12) E ( k x ,k z )= µ , (3)(see Supplementary Note 1 for details). We note that theFermi surface depends on the surface termination, and sodoes the OEE. We also confirm the surface dependencefrom numerical calculations as shown in Figs. 2d-i. Thisformula applies to any topological insulators such as Z -TIs (see Supplementary Note 2). Surface theory of OEE for a cylinder.
From thisslab calculation, we calculate the OEE for a cylinder ge-ometry. We consider a current along the z direction in aone-dimensional quadrangular prism with xz and yz sur-faces (surfaces I-IV in Fig. 3) through its surface Hamil-tonian. Let L x and L y denote the system sizes along the x and y directions, respectively. Because the OEE is sen-sitive to differences in crystal surfaces, as shown in slabsystems, we consider the individual surfaces separately.In particular, in Chern insulators we can calculate theenergy eigenstates for the whole system from those for the a b y x0 L y L x j x (cid:127) =v x (cid:127) |u xz (cid:127) | j y (cid:127)(cid:127) =v y (cid:127)(cid:127) |u yz (cid:127) | xyz (cid:127)V (cid:127)(cid:127)(cid:127)(cid:127)(cid:127)(cid:127) y (cid:127)V L x L y (cid:127)(cid:127)(cid:127) FIG. 3:
One-dimensional prism. ( a ) Cross section ofa one-dimensional system along xy -plane. ( b ) Currentconservation at the corners.surface Hamiltonians. For simplicity, we assume twofoldrotation symmetry C z of the system, which relates be-tween I and III, and between II and IV. Then, only thesurface I and II are independent. We write down theeigenequations for these surfaces as H I ψ I k x k z ( x, z ) = E I k x k z ψ I k x k z ( x, z ) , (4) H II ψ II k y k z ( y, z ) = E II k y k z ψ II k y k z ( y, z ) , (5)where H I and H II are the surface Hamiltoni-ans for the surfaces I and II, respectively, and ψ I k x k z = u I k x k z ( k x , k z ) e ik x x e ik z z and ψ II k y k z = u II k y k z ( k y , k z ) e ik y y e ik z z are Bloch eigenstates on thesurfaces I and II, respectively. We can determine theseeigenstates from four conditions, equality of the energyeigenvalues, current conservation at the corner ,periodic boundary condition on the crystal surface andthe normalization condition (see Supplementary Note3).Thus we obtain a formula for OEE in a one-dimensional prism of a three-dimensional Chern insulator M OEE z = − e τ E z ¯ h π ) Z π/c − π/c dk z L x ∂k I x ∂k z + L y ∂k II y ∂k z L x v I x + L y v II y (cid:12)(cid:12)(cid:12) E = µ , (6)where v I x = h ∂E I ∂k x and v II y = h ∂E II ∂k y and k I x ( k z , E ) and k II y ( k z , E ) are functions obtained from E = E I k x k z and E = E II k y k z , respectively. When v x and v y are almostindependent of k z , we approximate equation (6): M OEE z = 2 L x h v I x i M I , OEE z, slab + L y h v II y i M II , OEE z, slab L x h v I x i + L y h v II y i , (7)where M I , OEE z, slab and M II , OEE z, slab represent the OEE for a slab(equation (3)) with the surface I and that with the surfaceII, respectively. Thus, the OEE of the one-dimensionalsystem can be well approximated by equation (7) ex-pressed in terms of that for the slabs along xz and along yz planes. In general topological insulators, we can also deriveOEE in terms of a simple picture of a combined circuit,consisting of four surfaces I- IV with anisotropic transportcoefficients. We obtain M OEE z = j circ = L x σ I xx σ I xz + L y σ II yy σ II yzL x σ I xx + L y σ II yy E z , (8)where j circ is the circulating current density withinthe xy plane around the prism per unit length alongthe z -direction. σ I , II ij is the electric conductivity ten-sor for the surfaces I and II (see Supplementary Note4). On the other hand, we can also show M I , OEE z, slab = σ I xz E z , M II , OEE z, slab = σ II yz E z . In Chern insulators, by us-ing σ I xx ∝ h v x i and σ II yy ∝ h v y i , we arrive at equation(7). Thus, we can calculate the OEE from the surfaceelectrical conductivity from equation (8), which dependson the aspect ratio L x /L y .We numerically comfirm that the results of direct cal-culation by equation (2) and those for surface calculationby equation (7) agree well (Fig. 4a-c). When the inter-layer hopping is large (Fig. 4c), they slightly deviate fromeach other. This is because we cannot ignore the k z de-pendence of v x and v y and they are out of the scope ofthe approximate expression (7). Finite-size effect.
In our approximation theory, we as-sumed that the surface current is localized at the out-ermost sites and ignored a finite penetration depth. Infact, we can fit well the data with various system sizeswith a trial fitting function which includes a finite-sizeeffect in equation (7) (see Methods and SupplementaryNote 5). From these results, the finite-size effect is of theorder 1 /L in the leading order, coming from the finitepenetration depth. When the system size is much largerthan the penetration depth, the result is well describedby the surface theory as shown in Fig. 4d. Material.
Topological insulators without inversion sym-metry can be a good platform for obtaining large OEE,because the current flows on the surface. Therefore, theclosed loop created by the current is macroscopic and itefficiently induces the orbital magnetization, in contrastto the conventional OEE in metals, where sizes of cur-rent loops are microscopic. Moreover, the surface statesof topological materials are robust against perturbationscaused by impurities.Under the non-inversion-symmetry constraint, we can-not diagnose Z -TIs easily because the Z topologicalinvariant is expressed in terms of k -space integrals. Ouridea here is to use S symmetry to diagnose Z -TIs, wherewe only need to calculate wavefunctions at four momentaaccording to the symmetry-based indicator theories .After searching in the topological material database , wenotice that Cu ZnSnSe with and CdGeAs with are two ideal candidates of Z -TIs with a direct gapfor obtaining a large OEE (see Supplementary Note 6 fordetails). In the following, we will use Cu ZnSnSe , whichonly has S symmetry as shown in Fig. 5a, as an exampleto show the magnitude of the OEE with different surfacesand different surface terminations. a bdc Surface HamiltonianTight-binding modelL x /L y Fitting L x =10a, L y =30aL x =20a, L y =30aL x =30a, L y =30a μ/t y L x =10a, L y =30aL x =20a, L y =30aL x =30a, L y =30a μ/t y L x =10a, L y =30aL x =20a, L y =30aL x =30a, L y =30a μ/t y FIG. 4:
OEE in one-dimensional systems. ( a-c )OEE calculated from two different methods; one is adirect calculation by equation (2) (solid lines) and theother is by a combination of calculation results forsurfaces along xz and yz planes based on equation (7)(dashed lines). Parameter values are ( a ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y ,( b ) t x = b x = 1 . t y , t y = b y = m , t = 0 . t y and t = 0 . t y and ( c ) t x = t y = m = b x = b y , t = 0 . t y and t = 0 . t y . ( d ) Dependence on the aspect ratiowithin the xy plane with parameters t x = b x = 1 . t y , t y = b y = m , t = 0 . t y , t = 0 . t y and µ = 0. Bluepoints represent the result of equation (2) for varioussystem sizes. With parameter values t x = b x = 1 . t y , t y = b y = m , t = 0 . t y , t = 0 . t y and µ = 0, wecalculate the OEE with various system sizes fromequation (2). The system sizes are( L x , L y ) ∈ { a, a, . . . , a } × { a, a, . . . , a } .Red points represent the result of fitting the numericalresults of equation (2) with the fitting function equation(13), and its limit for L x , L y → ∞ is shown as thedashed line. The solid line represent the results ofequation (7).Since the magnetoelectric tensor for the space group , defined by M = αE , is α = α α α − α
00 0 0 , wecan obtain an orbital magnetization M OEE1 by adding anexternal electric field E , through the surface currentsboth on the [001] surface and on the [010] surface thanksto the nonzero α (see Supplementary Note 7 for de-tails).Figure 5c is the band structure of Cu ZnSnSe withgap, through the first-principle calculations whose detailsare explained in Methods. On the [001] surface, termi-nations with Cu-Sn layer (surface A) and with Se layer(surface B) have different surface energies and Fermi sur-faces, as shown in Fig. 5d-g, which contribute to a mag-netoelectric susceptibility of α A11 = − . × s − Ω − · τ and α B11 = − . × s − Ω − · τ , respectively. Onthe A surface, there is a single surface Dirac cone at Γpoint, forming an electron-like Fermi surface. On theB surface, the Dirac cone at Γ point forms an almostzero Fermi surface, but two surface Dirac cones at two¯ X momenta form two hole-like Femi surfaces. Becausethe Fermi surfaces on the B surface are much larger thanthose on the A surface, the magnetoelectric susceptibilityon the B surface is one order of magnitude larger thanthat on the A surface. Similar calculations on the [010]surface are in the Supplementary Note 8, and the resultis α C11 = − . × s − Ω − · τ for the surface C.Let us compare the results with metallic materials inthe bulk. For simplicity, we focus on the cases with theelectric field E and the resulting magnetization M OEE along the z direction. This kinetic magnetoeletric re-sponse is expressed as M OEE z = α zz E z , and the conduc-tivity is j z = σ zz E z . Thus the magnetization in responseto the current is M OEE z = ( α zz /σ zz ) j z . In the relaxationtime approximation, both α zz and σ zz are proportionalto the relaxatoin time τ . For simplicity we consider thesystem to be a cube with its size L × L × L . In the bulkmetallic systems, the current is carried by the bulk states,and α zz ∝ L , σ zz ∝ L . On the other hand, in topo-logical systems, only the surface conducts the current,and the conductivity σ zz scales as σ zz ∝ L − . On theother hand, we have shown α zz ∝ L , which means that α zz is an intensive quantity. Thus, in topological insula-tors, the scaling of the OEE as a response to the electricfield is represented by the response coefficient α zz ∝ L .Meanwhile, the response coefficient α zz /σ zz to the cur-rent is proportional to L . It means that as a responseto the current, topological materials will generate a largeamount of orbital magnetization as compared to metals.We compare our results with OEE in p-doped tel-lurium, which has chiral crystal structure . For the ac-ceptor concentration N a = 4 · cm − at 50K, theinduced orbital magnetizations is M OEE z = 7 . · − µ B / atom ∼ . × − A/m by a current density j z =1000A/cm . Thus the response coefficient of the or-bital magnetization M OEE z to the current density j z is α zz /σ zz = 1 . × − m. If we approximate σ zz by σ zz ∼ N a e τm with the electronic charge e and mass m ,we get α = 2 . × s − Ω − · τ . For other acceptor con-centrations N a = 4 · cm − and N a = 1 · cm − ,one can similarly get α zz = 2 . × s − Ω − · τ and α zz = 2 . × s − Ω − · τ . Thus, the size of α zz for thetopological insulator Cu ZnSnSe is larger than that ofTe by two to five orders of magnitude.On the other hand, in topological insulators, the in-duced orbital magnetization as a response to the currentbecomes huge compared with metals. To show this, weconsider a system with surfaces having anisotropic trans-port coefficients. σ xz /σ zz = tan θ , where θ describesan angle between the electric field along the z directionand the surface current density j surf . For example, forthe [001] surface of Cu ZnSnSe , we get σ A21 = 2 α A11 = − . × s − Ω − · τ , σ A11 = 7 . × s − Ω − · τ , σ B21 =2 α B11 = − . × s − Ω − · τ , σ B11 = 2 . × s − Ω − · τ ,which yield tan θ A = − .
49 and tan θ B = − .
26 by iden-tifying x = 2 and z = 1. Then the total current alongthe z direction is 4 Lj z while the circulating current is j circ = j surf x = j surf z tan θ . Thus the magnetization re-sponse M OEE z to the current density j z (= 4 Lj surf z /L ) is M OEE z /j = ( j circ /j surf z )( L/
4) = ( L/
4) tan θ . Thus for themacroscopic system size, the response M OEE z /j z is alsoof the macroscopic size, and it is many order of mag-nitude larger than that in tellurium, where M OEE z /j z isevaluated to be M OEE z /j z = 1 . × − m. Discussion
In summary, we propose OEE in topological insulatorswith chiral structure. This OEE is carried by surfacecurrent due to the asymmetric crystal structure of thesurface. Therefore, the OEE is sensitive to surface ter-minations, and it cannot be defined as a bulk quantity.We derive a formula for the OEE as a surface quantityusing the surface Hamiltonian, and show that it fits withnumerical results.In theoretical treatments, atomic orbitals can classifythe orbital magnetization into intraatom and interatomcontributions. Some atomic orbitals such as p x ± ip y haveorbital angular momentum, which leads to correspondingintraatomic orbital magnetization. On the other hand,the hopping between atoms lead, to the interatomic or-bital magnetization. In tight-binding models with atomicorbitals, they are separately calculated. For example, inRefs. 10 and 11, the intraatomic orbital magnetization isstudied, while Refs. 1, 2, and 15 consider the interatomiccontributions. In real mateirals, these two contributionsare not separable, and in the ab initio calculation , theirsum is calculated. In this paper, we found that in topo-logical materials, the interatomic contribution is muchlarger due to the macroscopic current loop. We showthat the response to the current in topological insulatorsis much larger than in matels. MethodsDetails of the first-principle calculations
First-principle calculations of Cu ZnSnSe are implementedin the Vienna ab initio simulation package (VASP) with Perdew-Burke-Ernzerhof exchange correlation. A Γ-centered Monkhorst-Pack grid with 10 × × k -pointsand 460.8 eV for the cut-off energy of the plane wave ba-sis set is used for the self-consistent calculation. Surfacestates and Fermi surfaces calculations are performed bythe tight-binding model obtained by the maximally lo-calized Wannier functions . Details of the model Hamiltonian.
We consider aChern insulating system with a chiral crystal structure.The model is composed of infinite layers of the two-dimensional Wilson-Dirac model . The lattice sitesare expressed by ( i, j, l ), with i , j , l being integers, spec-ifying the x , y and z -coordinates. At each lattice site,we consider two orbitals 1 and 2. Let c i,j,l,σ denote theannihilation operator of electrons at the ( i, j, l )-site withorbital σ (= 1 , c i,j,l = ( c i,j,l, , c i,j,l, ) T .The model Hamiltonian is H = H WD ( m, t x , t y , b x , b y ) + H interlayer ( t , t ) , where H WD is an in-plane Wilson-Dirac Hamiltonian, and H interlayer is an interlayer Hamil- tonian representing a structure similar to right handedsolenoids. The in-plane Wilson-Dirac Hamiltonian is H WD = m X i,j,l c † i,j,l σ z c i,j,l − it x X i,j,l ( c † i,j,l σ x c i +1 ,j,l − H.c. ) − it y X i,j,l ( c † i,j,l σ y c i,j +1 ,l − H.c. )+ b x X i,j,l ( c † i,j,l σ z c i +1 ,j,l + H.c. − c † i,j,l σ z c i,j,l )+ b y X i,j,l ( c † i,j,l σ z c i,j +1 ,l + H.c. − c † i,j,l σ z c i,j,l ) , (9)where H.c. stands for Hermitian conjugate of the pre-ceding terms, † represents Hermitian conjugate, and m, t x , t y , b x and b y are real parameters. This Hamilto-nian H WD can be rewritten in the momentum space as˜ H WD ( k ) = t x sin k x aσ x + t y sin k y aσ y + ( m − b x (1 − cos k x a ) − b y (1 − cos k y a )) σ z , (10)where k is the Bloch wavenumber. An isotropic ver-sion of the two-dimensional Wilson-Dirac model with b ≡ b x = b y and t x = t y exhibits the Chern insulat-ing phase when 0 < m/b < < m/b < . Nextwe add interlayer hoppings, including a direct hopping t along the z -axis and a chiral hopping t , where t and t are real parameters. To describe the chiral hopping t , the lattice sites in the square lattice in each layer intogroups of four sites, (2 i − , j − , (2 i − , j ) , (2 i, j − i, j ) where i and j are integers, and we introducechiral hoppings between the groups on the neighboringlayers. Then the total Hamiltonian for this model on atetragonal lattice is given by H = H WD + H interlayer , (11)where H interlayer = t X i,j,l ( c † i − , j − ,l c i, j − ,l +1 + H.c. )+ t X i,j,l ( c † i, j − ,l c i, j,l +1 + H.c. )+ t X i,j,l ( c † i, j,l c i − , j,l +1 + H.c. )+ t X i,j,l ( c † i − , j,l c i − , j − ,l +1 + H.c. )+ t X i,j,l ( c † i,j,l c i,j,l +1 + H.c. ) . (12)These hoppings in H interlayer form structures similar toright-handed solenoids. When H WD in the Chern insu-lator phase, even if H WD is perturbed by H interlayer , thesystem remains in the Chern insulator with the Chern SnZnCuSe b P (a)(b) ГM Nb Xb Г X N Г M(c)210-1-2 E ne r g y ( e V ) ГX –– X–M– -0.4-0.3-0.2-0.10.20.10.0 (-",-") k k (- " , " ) (",-") X–Г– M–[001] (-0.03",-0.03") k k (d) (e)(f) (g) (0.03",-0.03") (- . " , . " ) -0.4-0.3-0.2-0.10.20.10.0 ( . " , . " )( " , " ) [001] surface A[001] surface B Fermi arc on surface AFermi arc on surface B FIG. 5:
First principle calculations on Cu ZnSnSe . ( a ) Crystal structure of Cu ZnSnSe . ( b ) Brillouin zoneand surface Brillouin zone along [001] direction. ( c ) Electronic structure with spin-orbit coupling for the bulk. ( d-e )Surface states and Fermi surface calculation on the [001] surface with Cu-Sn layer termination (surface A). ( f-g )Surface states and Fermi arcs calculation on the [001] surface with Se layer termination (surface B).number within the xy plane equal to − t and t are small. In the main text, we are interestedin the OEE due to the topological surface states in thetopological Chern insulating phase, in which the Fermienergy is in the energy gap. Fitting function for OEE.
By taking into account thefinite-size effect in equation (7), we give a fitting fuction.The finite penetration depth of the surface states willlead to O (1 /L ) correction to the OEE, and that aroundthe corner will lead to O (1 /L ) correction . Thus, thefitting function is M OEE z = w L x + w L y + w + w L x + w L y w L x + w L y , (13)where w i ( i = 1 , , . . . ,
7) are real constants.
Acknowledgement
This work was supported by Japan Society for thePromotion of Science (JSPS) KAKENHI Grants No.JP18H03678, and No. JP20H04633, and by Elements strategy Initiative to Form Core Research Center (TIES),from MEXT Grant Number JP-MXP0112101001.
Additional information
The authors declare no competing financial interests.
Data availability statement
The datasets generated during and/or analysed duringthe current study are available from the correspondingauthors on reasonable request.
Code availability statement
The source code for the calculations performed in thiswork is available from the corresponding authors uponreasonable request.
Author contribution
All authors contributed to the main contents of this work.K. O. performed the model calculation and formulatedthe theory for topological systems, through the discus-sions with T. Z. and S. M. T. Z. performed the ab initio calculation. S.M. conceived and supervised the project.All authors drafted the manuscript. ∗ [email protected] T. Yoda, T. Yokoyama, and S. Murakami, Sci. Rep ,12024 (2015). T. Yoda, T. Yokoyama, and S. Murakami, Nano Lett. ,916 (2018). T. Furukawa, Y. Shimokawa, K. Kobayashi, and T. Itou,Nat. Commun. , 954 (2017). S. Zhong, J. E. Moore, and I. Souza, Phys. Rev. Lett. ,077201 (2016). S. S. Tsirkin, P. A. Puente, and I. Souza, Phys. Rev. B , 035158 (2018). Y.-Q. Wang, T. Morimoto, and J. E. Moore, Phys. Rev.B , 174419 (2020). J. E. Moore and J. Orenstein, Phys. Rev. Lett. , 026805 (2010). I. Sodemann and L. Fu, Phys. Rev. Lett. , 216806(2015). Q. Ma, S.-Y. Xu, H. Shen, D. MacNeill, V. Fatemi, T.-R. Chang, A. M. M. Valdivia, S. Wu, Z. Du, C.-H. Hsu,S. Fang, Q. D. Gibson, K. Watanabe, T. Taniguchi, R. J.Cava, E. Kaxiras, H.-Z. Lu, H. Lin, L. Fu, N. Gedik, andP. Jarillo-Herrero, Nature , 337 (2018). V. A. Shalygin, A. N. Sofronov, L. E. Vorob’ev, and I. I.Farbshtein, Physics of the Solid State , 2362 (2012). T. Koretsune, R. Arita, and H. Aoki, Phys. Rev. B ,125207 (2012). T. Furukawa, Y. Watanabe, N. Ogasawara, K. Kobayashi,and T. Itou, “Chirality-induced electrical generation ofmagnetism in nonmagnetic elemental tellurium,” (2020),arXiv:2010.09210 [cond-mat.mtrl-sci]. J. Rou, C. S¸ahin, J. Ma, and D. A. Pesin, Phys. Rev. B , 035120 (2017). C. S¸ahin, J. Rou, J. Ma, and D. A. Pesin, Phys. Rev. B , 205206 (2018). D. Hara, M. S. Bahramy, and S. Murakami, Phys. Rev. B , 184404 (2020). W.-Y. He, D. Goldhaber-Gordon, and K. T. Law, Nat.Commun. , 1650 (2020). V. Edelstein, Solid State Commun. , 233 (1990). E. L. Ivchenko and G. E. Pikus, ZhETF Pisma Redaktsiiu , 640 (1978). L. Levitov, Sov. Phys. JETP , 133 (1985). R. Naaman and D. H. Waldeck, The Journal of PhysicalChemistry Letters , 2178 (2012). B. Gohler, V. Hamelbeck, T. Z. Markus, M. Kettner, G. F.Hanne, Z. Vager, R. Naaman, and H. Zacharias, Science , 894 (2011). M. Kettner, V. V. Maslyuk, D. N¨urenberg, J. Seibel,R. Gutierrez, G. Cuniberti, K.-H. Ernst, and H. Zacharias,The Journal of Physical Chemistry Letters , 2025 (2018). B. P. Bloom, V. Kiran, V. Varade, R. Naaman, and D. H.Waldeck, Nano Letters , 4583 (2016). O. B. Dor, N. Morali, S. Yochelis, L. T. Baczewski, andY. Paltiel, Nano Letters , 6042 (2014). T. Kimura, T. Goto, H. Shintani, K. Ishizaka, T. Arima,and Y. Tokura, Nature , 55 (2003). H. Katsura, N. Nagaosa, and A. V. Balatsky, Phys. Rev.Lett. , 057205 (2005). M. Fiebig, Journal of Physics D: Applied Physics , R123(2005). N. A. Spaldin and M. Fiebig, Science , 391 (2005). F. Wilczek, Physical Review Letters , 1799 (1987). A. M. Essin, J. E. Moore, and D. Vanderbilt, Phys. Rev. Lett. , 146805 (2009). A. M. Essin, A. M. Turner, J. E. Moore, and D. Vander-bilt, Phys. Rev. B , 205104 (2010). A. Malashevich, I. Souza, S. Coh, and D. Vanderbilt, NewJ. Phys. , 053032 (2010). S. Coh, D. Vanderbilt, A. Malashevich, and I. Souza, Phys.Rev. B , 085108 (2011). D. Ceresoli, T. Thonhauser, D. Vanderbilt, and R. Resta,Phys. Rev. B , 024408 (2006). T. Thonhauser, D. Ceresoli, D. Vanderbilt, and R. Resta,Phys. Rev. Lett. , 137205 (2005). D. Xiao, J. Shi, and Q. Niu, Phys. Rev. Lett. , 137204(2005). K.-T. Chen and P. A. Lee, Phys. Rev. B , 195111 (2012). R. Bianco and R. Resta, Phys. Rev. Lett. , 087202(2013). A. Marrazzo and R. Resta, Phys. Rev. Lett. , 137201(2016). L. Fu, C. L. Kane, and E. J. Mele, Phys. Rev. Lett. ,106803 (2007). J. E. Moore and L. Balents, Phys. Rev. B , 121306(R)(2007). D. Hsieh, D. Qian, L. Wray, Y. Xia, Y. S. Hor, R. J. Cava,and M. Z. Hasan, Nature , 970 (2008). A. Raoux, M. Polini, R. Asgari, A. R. Hamilton, R. Fazio,and A. H. MacDonald, Phys. Rev. B , 073407 (2010). A. Concha and Z. Teˇsanovi´c, Phys. Rev. B , 033413(2010). R. Takahashi and S. Murakami, Phys. Rev. Lett. ,166805 (2011). H. C. Po, A. Vishwanath, and H. Watanabe, Nature com-munications , 1 (2017). Z. Song, T. Zhang, Z. Fang, and C. Fang, Nature commu-nications , 1 (2018). E. Khalaf, H. C. Po, A. Vishwanath, and H. Watanabe,Physical Review X , 031070 (2018). T. Zhang, Y. Jiang, Z. Song, H. Huang, Y. He, Z. Fang,H. Weng, and C. Fang, Nature , 475 (2019). L. Guen, W. Glaunsinger, and A. Wold, Mater. Res. Bull. , 463 (1979). G. Kresse and J. Hafner, Phys. Rev. B , 558 (1993). G. Kresse and J. Hafner, Phys. Rev. B , 14251 (1994). G. Kresse and J. Furthmuller, Comput. Mater. Sci. , 15(1996). A. A. Mostofi, J. R. Yates, Y.-S. Lee, I. Souza, D. Vander-bilt, and N. Marzari, Comput. Phys. Commun. , 685(2008). M. Creutz and I. Horv´ath, Phys Rev D , 2297 (1994). Y. Yoshimura, K.-I. Imura, T. Fukui, and Y. Hatsugai,Phys. Rev. B , 155443 (2014). rbital Edelsein e ff ect in topological materials: Supplemental material Ken Osumi, Tiantian Zhang,
1, 2 and Shuichi Murakami
1, 2, ∗ Department of Physics, Tokyo Institute of Technology, Tokyo 152-8551, Japan TIES, Tokyo Institute of Technology, Tokyo 152-8551, Japan (Dated: today) a r X i v : . [ c ond - m a t . m e s - h a ll ] F e b upplementary Note 1. Surface theory in a slab of a 3D Chern insulator In this section, for a slab of a 3D Chern insulator we explain how to rewrite Eq. (2) into Eq. (3)in the main text. Let us set the Chern number along the xy plane to be C xy = −
1. The energy ofthe chiral surface states on the two surfaces y = y ± based only on an intralayer Hamiltonian withinthe xy plane, representing the Chern insulator, is E ± = ∓ ~ v x k x ( v > , where k x is the wavevectoralong the x -direction and v x is the velocity along the x -direction. It means that the surface currentis along the x -direction. When we include the interlayer Hamiltonian having chiral hoppings, thesurface current on the surface will acquire a nonzero z -component as shown in Fig. 2a. Therefore,the energies of the chiral surface states of this model on the surface y = y ± are E ± ( k x , k z ) = ∓ ~ v x k x + ~ v z k z . (S1)Here, the velocity v z along the z -direction is added to represent chiral nature of the hopping alongthe surface. This simple form means that the velocity of the chiral surface states is a constantvalue ( v x , v z ). Nonetheless, this formula needs to be modified by the following reason. The energyshould be a periodic function of k z with a period of the reciprocal lattice vector G z = π c along the z -direction: E ( k x , k z ) = E ( k x , k z + π c ) , but Eq. (S1) does not satisfy this requirement and we cannotuse Eq. (S1). On the other hand, since E ( k x , k z ) represents the topological chiral surface states, itis not periodic along the k x direction. For example, on the y = y − surface, by increasing k x thesurface states come out of the bulk valence band, and finally go into the bulk conduction band, asshown in Fig. 2b. Let E − = E − ( k x , k z ) be the surface-state dispersion on the y = y − surface. Forsimplicity, we assume C z symmetry of the system. Then the surface state dispersion on the y = y + surface is given by E + = E + ( k x , k z ) = E − ( − k x , k z ) . (S2)Let ψ + ( − ) ( k x , k z ) denote a wave function on the y = y + ( y = y − ) surface of the crystal. In thiscalculation, we assume that the surface states are sharply localized at y = y ± , which means that weignore finite-size e ff ects due to a finite penetration depth. From Eq. (3), we evaluate the OEE as M OEE z , slab = M z , + + M z , − , (S3) M z , + = π ) Z π/ c − π/ c dk z L y Z π/ a − π/ a dk x e τ E z ~ ∂ f ( E + ( ˜k )) ∂ k z y + v x , + ( ˜k ) , (S4) M z , − = π ) Z π/ c − π/ c dk z L y Z π/ a − π/ a dk x e τ E z ~ ∂ f ( E − ( ˜k )) ∂ k z y − v x , − ( ˜k ) , (S5)2here ˜k = ( k x , k z ) and v x , ± ( ˜k ) = ~ ∂ E ± ( ˜k ) ∂ k x . Therefore, by noting that L y = y + − y − and v x , + ( ˜k ) = − v x , − ( ˜k ), we obtain the OEE as M OEE z , slab = − e τ E z ~ π ) Z BZ d ˜k ∂ f ( E − ( ˜k )) ∂ k z v x , − ( ˜k ) . (S6)Henceforth we omit the subscript “ − ” in E − . In the limit of zero temperature T → , the k z derivative of the Fermi-Dirac function is express as ∂ f ( E ( ˜k )) ∂ k z = − δ ( E − µ ) ∂ E ( ˜k ) ∂ k z , (S7)which allows us to rewrite Supplementary Eq. (S6). Thus, in terms of the surface states at theFermi energy, we can rewrite Supplementary Eq. (S6) as M OEE z , slab = e τ E z ~ π ) Z π/ c − π/ c dk z Z E u ( k z ) E l ( k z ) dE δ ( E − µ ) ∂ E ( k x , k z ) ∂ k z , (S8)where E l ( k z ) and E u ( k z ) are the lower and upper bounds of the surface state energies at k z . Byperforming the integral over E , we obtain M OEE z , slab = e τ E z ~ π ) Z π/ c − π/ c dk z ∂ E ( k x , k z ) ∂ k z (cid:12)(cid:12)(cid:12)(cid:12) E ( k x , k z ) = µ . (S9)This is Eq. (3) in the main text for Chern insulators, with sgn (cid:16) ∂ E ( k x , k z ) ∂ k x (cid:17) = + E ( k x , k z ) = µ (constant) yields ∂ E ∂ k z + dk x dk z (cid:12)(cid:12)(cid:12)(cid:12) E = µ ∂ E ∂ k x = , (S10)we can rewrite Eq. (S9) as M OEE z , slab = − e τ E z ~ π ) Z π/ c − π/ c dk z v x ( k FS x ( k z ) , k z ) dk FS x ( k z ) dk z , (S11)where k FS x is a function of k z defined by E ( k FS x ( k z ) , k z ) = µ . Supplementary Note 2. Extension of the theory to Z -TIs In Chern insulators, the topological surface states are chiral, and their energy can be amonotonous function of the wavenumber. In contrast, in Z -TIs, the eigenenergy of the topo-logical surface states is not monotonous. Nonetheless, by dividing the Fermi surface, we cancalculate the OEE in the same way as in the calculation of Chern insulators. For example, in the3 x k z Fermi surface k x1 (k z )k x2 (k z )k z+ k z- O FIG. S1.
Fermi surface on the y = y − surface in the STI. k x and k x are the value of k x ( k x > k x )on the Fermi surface. Inside and outside the Fermi surface represent the regions of E − ( k x , k z ) < µ and E − ( k x , k z ) > µ , respectively. case of a strong topological insulator (STI) with surface Fermi surfaces given in SupplementaryFig. 1, the contributions from k x ( k z ) and k x ( k z ) are written from Eq. (S8) as M OEE z , slab , k x = e τ E z ~ π ) Z k + z k − z dk z Z E b ( k z ) E a ( k z ) dE δ ( E − µ ) ∂ E ( k x , k z ) ∂ k z (cid:12)(cid:12)(cid:12)(cid:12) k x = k x ( k z ) = e τ E z ~ Z π/ c − π/ c dk z (2 π ) ∂ E ( k x , k z ) ∂ k z (cid:12)(cid:12)(cid:12)(cid:12) k x = k x ( k z ) , (S12) M OEE z , slab , k x = e τ E z ~ π ) Z k + z k − z dk z Z E a ( k z ) E b ( k z ) dE δ ( E − µ ) ∂ E ( k x , k z ) ∂ k z (cid:12)(cid:12)(cid:12)(cid:12) k x = k x ( k z ) = − e τ E z ~ Z π/ c − π/ c dk z (2 π ) ∂ E ( k x , k z ) ∂ k z (cid:12)(cid:12)(cid:12)(cid:12) k x = k x ( k z ) , (S13)where E ( k x ( k z ) , k z ) = µ, (S14) E ( k x ( k z ) , k z ) = µ, (S15) E a ( k z ) < µ < E b ( k z ) . (S16)The sum of Eqs. (S12) and (S13) is the OEE to STIs. Therefore, we obtain the OEE in Z -TIs asthe sum of these results: M OEE z , slab = e τ E z ~ Z π/ c − π/ c dk z (2 π ) ∂ E ( k x , k z ) ∂ k z sgn ∂ E ( k x , k z ) ∂ k x ! (cid:12)(cid:12)(cid:12)(cid:12) E ( k x , k z ) = µ . α zz , slab for the slab, defined by M OEE z , slab = α zz , slab E z , is givenby α OEE zz , slab = e τ ~ Z π/ c − π/ c dk z (2 π ) ∂ E ( k x , k z ) ∂ k z sgn ∂ E ( k x , k z ) ∂ k x ! (cid:12)(cid:12)(cid:12)(cid:12) E ( k x , k z ) = µ . (S17)From the above derivation, we can also directly derive that α zz , slab is related to the o ff -diagonalconductivity σ xz on the y = y − surface via σ zx = σ xz = α zz , slab = e τ ~ Z π/ c − π/ c dk z (2 π ) ∂ E ( k x , k z ) ∂ k z sgn ∂ E ( k x , k z ) ∂ k x ! (cid:12)(cid:12)(cid:12)(cid:12) E ( k x , k z ) = µ . (S18)Likewise, the longitudinal conductivity is written as σ zz = e τ ~ Z π/ c − π/ c dk z (2 π ) ∂ E ( k x , k z ) ∂ k z sgn ∂ E ( k x , k z ) ∂ k z ! (cid:12)(cid:12)(cid:12)(cid:12) E ( k x , k z ) = µ . (S19) Supplementary Note 3. Surface Theory in a one-dimensional prism of 3D Chern insulator
In this section, for a one-dimensional prism we explain how to rewrite Eq. (2) into Eq. (6). Wecan determine the eigenstates on the surfaces I and II, u I k x k z and u II k x k z from four conditions. The firstcondition is that the energy eigenvalues are equal: E I k x k z = E II k y k z . The second one is the currentconservation at the corner [1–3]: v I x ( k x , k z ) | u I k x k z | = v II y ( k y , k z ) | u II k y k z | , (S20)where v I x ( k x , k z ) = ~ ∂ E I k x k z ∂ k x , (S21) v II y ( k y , k z ) = ~ ∂ E II k y k z ∂ k y . (S22)The third one is the periodic boundary condition on the crystal surface:2 k x L x + k y L y = π n ( n = , ± , ± , . . . ) , (S23)where L x and L y are the lengths of the crystal in the x and y directions. This integer n canbe regarded as a quantum number labeling the eigenstates, and we write their eigenvalues as E ( n , k z )( = E I k x k z = E II k y k z ). The last one is the normalization condition:2 L x | u I k x k z | + L y | u II k y k z | = . (S24)5rom these conditions, u I k x k z and u II k y k z are | u I k x k z | = v I x L x v I x + L y v II y , | u II k y k z | = v II y L x v I x + L y v II y . (S25)Then we obtain OEE from Eq. (3) and Eq. (S25): M OEE z = π Z π/ c − π/ c dk z X n e τ E z ~ ∂ f ( E ( n , k z )) ∂ k z (cid:18) − e (cid:19) L x v I x + L y v II y . (S26)When L x and L y are large, k x , k y and E ( n , k z ) are regarded as continuous variables. By rewritingthe summation over n to an integral over E and by using Eq. (S23), we obtain the OEE as M OEE z = − e τ E z ~ π Z π/ c − π/ c Z dEdk z ∂ f ( E ( n , k z )) ∂ k z dndE L x v I x + L y v II y = e τ E z ~ π ) Z π/ c − π/ c dk z ∂ E ( n , k z ) ∂ k z (cid:12)(cid:12)(cid:12) E ( n , k x ) = µ , (S27)where we used dndE = ~ π L x v I x + L y v II y ! , (S28) ∂ f ( E ( n , k z )) ∂ k z = − δ ( E − µ ) ∂ E ( n , k z ) ∂ k z (S29)From Eq. (S23), let k I0 x ( k z , E ) and k II0 y ( k z , E ) be the values of k x and k y for a given value of k z and E . Therefore they satisfy n = π ( L x k I0 x ( k z , E ) + L y k II0 y ( k z , E )) (S30)Therefore, we finally get: ∂ E ( n , k z ) ∂ k z + ∂ E ( n , k z ) ∂ n π L x ∂ k I0 x ( k z , E ) ∂ k z + L y ∂ k II0 y ( k z , E ) ∂ k z = , (S31)Therefore, we rewrite Eq. (S27) as M OEE z = − e τ E z ~ π ) Z π/ c − π/ c dk z L x ∂ k I x ∂ k z + L y ∂ k II y ∂ k z L x v I x + L y v II y (cid:12)(cid:12)(cid:12)(cid:12) E = µ , (S32)where v I x = ~ ∂ E I ∂ k x and v II y = ~ ∂ E II ∂ k y and k I x ( k z , E ) and k II y ( k z , E ) are functions obtained from E = E I k x k z and E = E II k y k z , respectively. This is the formula for OEE in a one-dimensional prism of a three-dimensional Chern insulator. Supplementary Note 4. OEE and electric conductivity tensors
6e can also derive OEE in a Chern insulating system by using electric conductivity tensors onsurfaces. On the surface I along the xz plane, the surface current density along xz plane, j x and j z in response to the electric field ( E x , E z ), is represented as j I x = σ I xx E I x + σ I xz E z , (S33) j I z = σ I zx E I x + σ I zz E z . (S34)where σ I i j is the electric conductivity tensor. Since we assumed that the system is invariant undertwofold rotation with respect to the z -axis, the surface III has the same tranport property with thesurface I, but with an opposite normal vector. Namely, it has j III x = σ I xx E III x − σ I xz E z , (S35) j III z = − σ I zx E III x + σ I zz E z . (S36)We note that the electric field along the z -axis, E z , is common among the four surfaces. Similarly,on the surface II along the yz plane, the surface current density along the yz plane, j y and j z inresponse the electric field ( E y , E z ), is written as j II y = σ II yy E II y + σ II yz E z , (S37) j II z = σ II zy E II y + σ II zz E z . (S38)Similarly, on the surface IV one has has j IV y = σ II yy E IV y − σ II yz E z , (S39) j IV z = − σ II zy E IV y + σ II zz E z . (S40)We note that the electric field satisfies I surface E · d l = L x ( E I x − E III x ) + L y ( E II y − E IV y ) = , (S41)and that the current around the crystal is conserved: j circ ≡ j I x = j II y = − j III x = − j IV y . (S42)Combining Eqs. (S33)-(S42), we obtain M OEE z = j circ = L x σ I xx σ I xz + L y σ II yy σ II yzL x σ I xx + L y σ II yy E z . upplementary Note 5. Finite-size e ff ect In this section, we introduce our fitting result in detail. We calculate the OEE with varioussystem sizes from Eq. (2). By using these results, we determine the parameters in the fittingfunction (see Methods). With the hopping parameters t x = b x = . t y , t y = b y = m , t = . t y , t = . t y and µ =
0, we get the results fitting as M OEE z K = . L x + . L y − . − . L x − − . L y L x + . L y × − , (S43)where K ≡ e τ E z t y / ~ , if we take the limit of L x → ∞ and L y → ∞ with the ratio L x / L y keptconstant, we get: M OEE z K = . L x + . L y L x + . L y × − . (S44)On the other hand, the result of Eq. (7) based on the surface theory is fitted as: M OEE z K = . L x + . L y L x + . L y × − . (S45)Apart from the coe ffi cients of L x in the numerators, Eqs. (S44) and (S45) agree well. Since thee ff ect of L x in the numerators of Egs. (S44) and (S45) is much smaller than that of L y , it is di ffi cultto numerically obtain the coe ffi cient of L x in the numerator with good accuracy. Thus, we concludethat these results agree well as shown in Fig. 4d. Supplementary Note 6. Cu ZnGeTe After searching in the topological material database [4], we notice that Cu ZnSnSe [5] with and CdGeAs with are two ideal candidates of Z -TIs with a direct gap . AlthoughCu ZnGeTe [5] is also classified to be a Z -TI according to the topological material database,there are several electron and hole pockets at the Fermi energy. Supplementary Note 7. Magnetoelectric tensor for the space group In the main text, the kinetic magnetoelectric tensor for the space group ( S ), defined by M = α E , is α = α α α − α
00 0 0 . (S46)As an example, we calculate α in a one-dimensional quadrangular prism of Cu ZnSnSe fromEq. (8). From S symmetry in Cu ZnSnSe , we get σ I11 = σ II22 and σ I13 = − σ II23 . Hence, when thequadrangular prism preserves the S symmetry, i.e. when L x = L y , we get α = L , L , α become nonzero as opposed toEq. (S46), because the S symmetry is broken by the shape of the sample.Thus, in general systems, if a component of the kinetic magnetoelectric tensor α i j is nonzerofrom symmetry analysis, it means that the corresponding OEE appears regardless of the crystalshape, reflecting the symmetry of the crystal. On the other hand, when α i j =
0, the OEE does notoccur when the crystal shape preserve the point-group symmetry of the crystal, but otherwise theOEE may occur depending on the crystal shape.
Supplementary Note 8. [010] surface states of Cu ZnSnSe In the main text we calculated the [001] surface of Cu ZnSnSe . In this Supplementary note wecalculate the surface of [010] direction, we also calculated the surface states and magnetoelectricsusceptibility with two di ff erent terminations, i.e., surface C with Cu-Sn layer termination andsurface D with Se layer termination, as shown in Fig. S2. On the surface D, which is the Se layertermination, the Fermi arc is very small and buried in the bulk states, so we only calculate the α C for illustration. On the surface C, there is big hole-like pocket formed by the Dirac cone at ¯ X ,which contributes to large magnetoelectric susceptibility of α C11 = − . × τ s − Ω − .9 nZnCuSe b P(a) (b)ГM Nb Xb X– Г–M–[010]ГX –– Z–M– -0.2-0.10.20.10.0 (0,0) k k ( , ! ) (2!,0) (c) (d) ( ! , ! ) [010] surface C Fermi arc on surface C (0,0) k k (f) Fermi arc on surface D ( , ! ) ( ! , ! ) (2!,0) (e) [010] surface D-0.2-0.10.20.10.0 ГX –– Z–M– E ne r g y ( e V ) E ne r g y ( e V ) Z– FIG. S2.
First principle calculations of Cu ZnSnSe on the [010] surface. ( a-b ) Crystal structure andBrillouin zone for the bulk and surface along [010] direction. ( c-d ) Surface states and Fermi arcs calcula-tion on the [010] surface with Cu-Sn layer termination (surface C). ( f-g ) Surface states and Fermi surfacecalculation on the [010] surface with Se layer termination (surface D), which has a very small Fermi surfaceburied in the bulk states. Supplementary References ∗ [email protected]
1] A. Raoux, M. Polini, R. Asgari, A. R. Hamilton, R. Fazio, and A. H. MacDonald,
Phys. Rev. B ,073407 (2010).[2] A. Concha and Z. Teˆsanovi´c, Phys. Rev. B , 033413 (2010).[3] R. Takahashi and S. Murakami, Phys. Rev. Lett. , 166805 (2011).[4] T. Zhang, Y. Jiang, Z. Song, H. Huang, Y. He, Z. Fang, H. Weng, and C. Fang,
Nature , 475 (2019).[5] L. Guen, W. Glaunsinger, A. Wold,
Mater. Res. Bull. , 463-467 (1979)., 463-467 (1979).