Tunable dynamical magnetoelectric effect in antiferromagnetic topological insulator MnBi$_2$Te$_4$ films

More than forty years ago, axion was postulated as an elementary particle with a low mass and weak interaction in particle physics to solve the strong $\mathcal{CP}$ (charge conjugation and parity) puzzle. Axions are also considered as a possible component of dark matter of the universe. However, the existence of axions in nature has not been confirmed. Interestingly, axions arise as pseudoscalar fields derived from the Chern-Simons theory in condensed matter physics. In antiferromagnetic insulators, the axion field can become dynamical induced by spin-wave excitations and exhibits rich exotic phenomena, such as, the chiral magnetic effect, axionic polariton and so on. However, the study of the dynamical axion field is rare due to the lack of real materials. Recently, MnBi$_2$Te$_4$ was discovered to be an antiferromagnetic topological insulator with a quantized axion field protected by the inversion symmetry $\mathcal{P}$ and the magnetic-crystalline symmetry $\mathcal{S}$. Here, we studied MnBi$_2$Te$_4$ films in which both the $\mathcal{P}$ and $\mathcal{S}$ symmetries are spontaneously broken and found that the dynamical axion field and largely tunable dynamical magnetoelectric effects can be realized through tuning the thickness of MnBi$_2$Te$_4$ films, the temperature and the element substitution. Our results open a broad avenue to study axion dynamics in antiferromagnetic topological insulator MnBi$_2$Te$_4$ and related materials, and also is hopeful to promote the research of dark matter.

In antiferromagnetic (AFM) insulators, the axion field θ can be a non-quantized value, once both T and P are broken. The spin-wave excitations can induce fluctuations of the axion field [29][30][31], called as the dynamical axion field θ(r, t) (DAF) which has spatial and temporal dependence. θ(r, t) can lead to rich dynamical magnetoelectric (ME) effects, for example, the dynamical chiral magnetic effect [32][33][34][35], anomalous Hall effect [32], axionic polariton [29] and nonlinear electromagnetic effect (c) Schematic of the P, S-breaking mass m5 in MnBi2Te4 thin films. It expects to increase when the thickness of thin films decreases due to the finite-size effect. [36,37]. In principle, generic AFM insulators may exhibit the DAFs, such as Cr 2 O 3 [15,38], but the DAFs turn out to be too weak to be detected [31]. Large DAF requires the topologically nontrivial AFM insulators hosting nonzero spin Chern numbers [39]. Waals(vdW) layered material Mn 2 Bi 2 Te 5 and the superlattice (MnBi 2 Te 4 ) 2 (Bi 2 Te 3 ) 1 were predicted to host large DAFs [31,39], denoted as DAF insulators. However, though material Mn 2 Bi 2 Te 5 has been successfully grown in experiments, the sample quality needs to be further improved. The experimental synthesis of the superlattice (MnBi 2 Te 4 ) 2 (Bi 2 Te 3 ) 1 has not been reported. Therefore, the urgent challenge is to discover DAF insulators.
Recently, bulk MnBi 2 Te 4 was discovered to be an AFM topological insulator and it hosts a quantized axion field θ = π protected by the inversion symmetry P and the magnetic-crystalline symmetry S = T τ 1/2 , where τ 1/2 is the half translation operator [22,23,[40][41][42][43][44]. Interestingly, both P and S are spontaneously broken in MnBi 2 Te 4 films with even septuple (SL) layers. In this work, our results indicate that the fluctuations of AFM order can lead to a dynamical ME effect within a large tunable range in MnBi 2 Te 4 thin films, induced by the P and S-breaking mass m 5 .
As shown in Fig. 1a, the antiferromagnetic topological insulator MnBi 2 Te 4 [22] consists of SLs coupled to each other by a Van der Waals-type interaction. It features an A-type AFM magnetic ground state with the out-ofplane easy axis. Although the magnetic order explicitly breaks the time-reversal symmetry T , the bulk MnBi 2 Te 4 possesses a combined magnetic-crystalline symmetry S = T τ 1/2 , marked in Fig. 1a. In addition, the inversion sym-metry P is preserved for the bulk MnBi 2 Te 4 with the inversion center at the Mn layer in the middle of each SL. Interestingly, both the S symmetry and the inversion symmetry P are explicitly broken in even-SL MnBi 2 Te 4 thin films due to the finite size, schematically shown in Fig. 1b. As a result, the restriction for the quantization of the axion filed θ is removed, which will lead to a tunable dynamical ME effect as we will demonstrate in details.
Effective model analysis. The topological nature and the physical properties of MnBi 2 Te 4 are determined by the electronic structure near the Γ point because the minimum energy gap is at the Γ point. Therefore, we consider the effective k·p model with the four low-lying states indicates the spin. The state |P 1 + z stands for the bonding states of two Bi layers and |P 1 + z for the anti-bonding state of two Te layers. We start from the low-energy effective Hamiltonian of bulk MnBi 2 Te 4 which hosts the symmetries S, P and C 3z . The four-band k · p model up to the third order of k is written down as [22], where where m 4 is the mass term m 4 < 0 (m 4 > 0) for with (without) the band inversion in the band structure, and the Dirac matrices are Γ 0,1,..., where τ x,y,z and σ x,y,z denote the Pauli matrices for the orbital and spin, respectively. The third order term H 3 is given explicitly by The detailed parameters can be obtained by fitting with first-principles band structure and have been provided in Ref. [22], as will be adopted here. Further for MnBi 2 Te 4 thin films, both the S symmetry and inversion symmetry P are broken, but the combined symmetry PT is preserved. Based on the symmetry analysis, a mass term m 5 Γ 5 to the leading order adds to the effective k · p model (Eq. 1) for MnBi 2 Te 4 thin films. The mass term m 5 Γ 5 can significantly change the axion field θ from the quantized value (π). Actually, the mass term m 5 results from the interplay between the AFM order and the finite-size effect, and it is expected to decreases with increasing thickness of the films, and finally vanishes for the bulk limit, as illustrated in Fig. 1c. This is in contrast to the three-dimensional (3D) DAF insulator [29,31,39] with a non-vanishing m 5 term in the bulk limit. Based on the effective k · p model (Eq. 1), the energy gap at the Γ point of MnBi 2 Te 4 thin films is found to be E gap = 2 m 2 4 + m 2 5 . When the strength of spin-orbit coupling (SOC) is tuned from zero to the realistic value, the band inversion is induced along with the mass term m 4 changing from m 4 > 0 to m 4 < 0, and the energy gap E gap = 2|m 5 | reaches the minimum when the band inversion starts (m 4 = 0). Based on this consideration, we can quantitatively characterize the mass m 5 of MnBi 2 Te 4 thin films.
As an example, the fat band structures of 12-SL MnBi 2 Te 4 with different SOC are calculated and shown in Fig. 2a. When SOC is weak, the band structure indicates an energy gap E gap without the band inversion. As expect, when increasing SOC, E gap becomes smaller, and then reaches the minimum around SOC∼ 0.9, and then becomes bigger with the band inversion. The SOC dependence of E gap is further extracted and shown in Fig. 2b. In Fig. 2c, it explicitly shows the mass m 5 as a function of the thickness of N bi-SL MnBi 2 Te 4 thin films (red line), which indeed gradually decreases with increasing N . For the comparison, we also present the result of m 5 for Mn 2 Bi 2 Te 5 thin films (blue line). It is obvious that for large N towards the bulk limit, m 5 approaches zero (a finite value) for MnBi 2 Te 4 (Mn 2 Bi 2 Te 5 ), as expected. It is worth to notice that m 5 of MnBi 2 Te 4 films is tunable to zero within a larger range than that of Mn 2 Bi 2 Te 5 films, as shown in Fig. 2c. Furthermore, it is numerically found that m 5 exhibits an approximately linear behavior against 1/N , namely, the inverse of the film thickness for both MnBi 2 Te 4 and Mn 2 Bi 2 Te 5 films, as shown in Fig. 2d. Our results provide a practical guideline to tune the mass term m 5 through varying the thickness of MnBi 2 Te 4 films.
Dynamical magnetoelectric effect. The application of an electric field can induce a magnetization, and reversely, the application of a magnetic field can lead to a polarization, which is the typical ME effect. In 3D topological insulators, such as Bi 2 Se 3 [45], the ME response is isotropic, which can be described by the effective axion action S θ with the quantized axion field θ = π, protected by the time reversal symmetry T . In AFM topological insulator MnBi 2 Te 4 , both the inversion symmetry P and S symmetry constrain the axion field θ to the quantized value of π. However, for even-SL MnBi 2 Te 4 thin films, both P and T symmetries are broken and the mass term m 5 is induced, thus rendering unquantized θ. It should be pointed out that due to the finite-size effect in MnBi 2 Te 4 thin films, the ME response is generically anisotropic [46]. Moreover, since m 5 is related to the inherently fluctuating AFM order, the axion field will be time dependent ∂ t θ(t) = 0, leading to the dynamical ME effect. In the following, we will show that the dynamical ME effect can be realized and drastically tuned by adjusting the thickness of MnBi 2 Te 4 thin films, the temperature and element substitutions.
Here, we consider applying a magnetic field B in the plane of MnBi 2 Te 4 thin films (e.g. the x axis in Fig. 1a), which will induce a charge polarization P in the same direction. To describe such a ME response, we need to calculate the ME coefficient α defined as P = −αB and it can be obtained from the quantity γ ≡ (1/d) through α = γ(e 2 /2h) [17,46]. Here, η(z) is a dimensionless function defined as where the current correlation function Π xy (z, z 1 ) is given by the Kubo formula, Here, |nk is the normalized Bloch wavefunction of the nth electron subband, and j 3D (k, z) = (e/ )∂ k H 3D (k, z) represents the 3D in-plane current density operator, where H 3D (k, z) is obtained by discretizing the system into a tight-binding model along the z direction and fixing the periodic boundary condition in the x-y plane to preserve k x and k y as good quantum numbers.
In the bulk limit of 3D topological insulators, the quantity γ is related to the axion field by γ = θ/π. According to Ref. [29], the mass term m 5 could lead to a correction of θ to the linear order. Considering the inherent AFM fluctuations, which changes m 5 by δm 5 , θ should be divided into a static part θ 0 and a dynamical part δθ(t). As a result, γ also consists of a static part γ 0 and a dynamical part δγ(t) given by δγ = δm 5 /g with the coefficient 1/g = ∂γ/∂m 5 to the linear order of δm 5 . We can see that 1/g is an important quantity for a large dynamical ME effect.
In Fig. 3a, we present the quantity γ as a function of the thickness of N bi-SLs MnBi 2 Te 4 thin films (red line). It can be seen that with increasing N , γ gradually increases and approaches 1 at the bulk limit (θ = π) with N → ∞. The growing deviation from γ = 1 with decreasing N results from the induced m 5 term by the increasing finite-size effect. We also calculated γ of Mn 2 Bi 2 Te 5 films (yellow line) for the comparison. On the whole, the γ has similar features for both MnBi 2 Te 4 and Mn 2 Bi 2 Te 5 films due to the finite-size effect. The key difference is that γ gradually approaches 1 for MnBi 2 Te 4 films at the bulk limit but a finite value (< 1) for for Mn 2 Bi 2 Te 5 films, which indicates the bulk MnBi 2 Te 4 has m 5 = 0 (θ = π) and the bulk Mn 2 Bi 2 Te 5 has m 5 = 0 (θ ∼ π). The dependence of γ on the mass term m 5 for MnBi 2 Te 4 films with different thickness (N = 4, 6, 8, 16), shown in Fig. 3b. When N is large, γ quickly increases within the small m 5 region, and inversely when N is small, γ slowly increases. With further increasing m 5 , γ becomes to linearly decrease for all films, which indicates a large m 5 tends to suppress γ.
As for the dynamical part of γ, we calculated the quantity 1/g on different m 5 for N bi-SL MnBi 2 Te 4 films, shown in Fig. 3d. It is found that 1/g increases with decreasing m 5 , and the maximum value appears when m 5 → 0. In order to obtain a large quantity 1/g, it needs a large N with a small m 5 . We know that the evolution of m 5 depends on the thickness N of MnBi 2 Te 4 films, shown in Figs. 2c and 2d. Based on realistic values of m 5 induced by the finite-size effect, we then calculated the evolution of 1/g with varying N in Fig. 3c. We can see that 1/g has a large value for the ultra thin limit, and quickly decreases with increasing the film thickness (N ), and then crosses zero from positive to negative, and gradually saturates. The positive large 1/g for the small N indicates the increase of γ, and the saturation of 1/g for the large N indicates the linear suppression, which is consistent with Fig. 3b. Importantly, stress that we can further tune m 5 → 0 through raising the temperature to AFM Néel temperature to obtain a drastically enhanced large 1/g of MnBi 2 Te 4 films within a large tunable range for experiments, see the shaded region in Fig. 3d. At last, we also investigate the dependence of the dynamical ME effect on the mass term m 4 , with m 4 < 0 (m 4 > 0) indicating a topologically nontrivial (trivial) 3D electric structure. We choose three representative m 4 to show γ and 1/g for different thickness N of films in Fig. 4. In the topologically trivial case (m 4 = 50 meV), γ ∼ 0 (θ = 0). In the topologically nontrivial case (m 4 < 0), it can be seen that a larger |m 4 | exhibits a larger γ, which approaches 1 (θ = π) for N → ∞. According to Fig. 4b, 1/g in the topologically trivial case is quite small and obviously negligible compared to that in the topologically nontrivial case. Interestingly, 1/g can reach a maximum value through tuning |m 4 | and the thickness N , which also indicates that a proper band inversion (|m 4 |) can give the largest 1/g, not the strongest band inversion. Since m 4 is related to the band inversion, we can experimentally tune it through element substitutions in Mn(Bi 1−x Sb x ) 2 (Te 1−y Se y ) 4 family materials.
In conclusion, our work presents that AFM topological insulator MnBi 2 Te 4 family films, which have been suc-cessfully synthesized and well studied, provide a perfect material platform to realize the dynamical ME effects. More importantly, the dynamical ME effect is tunable and even drastically enhanced through engineering the finite-size effect, the temperature and the element substitution. The exotic axion dynamics, including dynamical chiral magnetic effect, axionic polariton and nonlinear electromagnetic effect, should be experimentally realized in the proposed real-material platform. Our discovery in condensed matter physics might also be helpful to detect the real axion particle and the dark matter in the universe [47].
Methods. We performed the first-principles calculations for bulk MnBi 2 Te 4 and Mn 2 Bi 2 Te 5 of different spinorbit coupling through employing the Vienna ab-initio simulation package(VASP) [48,49] and the generalized gradient approximation (GGA) with the Perdew-Burke-Ernzerhof(PBE) [50,51] type exchange-correlation potential was adopted, with the energy cutoff fixed to 420 eV. The k-point sampling grid of Brillouin zone in the self-consistent process was a Γ-centered Monkhorst-Pack k-point mesh of 8 × 8 × 4, and a total energy tolerance 10 −7 eV was adopted for self-consistent convergence. For MnBi 2 Te 4 and Mn 2 Bi 2 Te 5 films, we employed the maximally localized Wanneir functions (MLWF) [52,53] to construct tight-binding Hamiltonians from the firstprinciple calculations. The Mn−d, Te−p and Bi−p orbitals were initialized for MLWFs by Wannier90 [54].