Observation of the universal magnetoelectric effect in a 3D topological insulator

The electrodynamics of topological insulators (TIs) is described by modified Maxwell's equations, which contain additional terms that couple an electric field to a magnetization and a magnetic field to a polarization of the medium, such that the coupling coefficient is quantized in odd multiples of α/4π per surface. Here we report on the observation of this so-called topological magnetoelectric effect. We use monochromatic terahertz (THz) spectroscopy of TI structures equipped with a semitransparent gate to selectively address surface states. In high external magnetic fields, we observe a universal Faraday rotation angle equal to the fine structure constant α=e2/2hc (in SI units) when a linearly polarized THz radiation of a certain frequency passes through the two surfaces of a strained HgTe 3D TI. These experiments give insight into axion electrodynamics of TIs and may potentially be used for a metrological definition of the three basic physical constants.

M axwell's equations are in the foundation of modern optical and electrical technologies. To apply Maxwell's equations in conventional matter, it is necessary to specify constituent relations, describing the polarization P c (E) and magnetization M c (B) as a function of the applied electric and magnetic fields, respectively. Soon after the theoretical prediction [1][2][3] and experimental discovery of two-dimensional (2D) and three-dimensional (3D) topological insulators (TIs) 4,5 , it has been recognized that the constituent relations in this new phase of quantum matter contain additional cross-terms P t (B) and M t (E) when time-reversal symmetry is weakly broken 6 .
Here N is an integer and aE1/137 is the fine structure constant. The derivation of equation (1) is based on the topological field theory of time-reversal invariant insulators 6 . Its intriguing consequences are the universal Faraday rotation angle y F j j¼ a, when a linearly polarized electromagnetic radiation passes through the top and bottom topological surfaces 6,7 , and magnetic monopole images, induced by electrical charges in proximity to a topological surface 8 . However, experimental verification of these topological magnetoelectric effects (TMEs) has been lacking. As the modified Maxwell's equations describing electrodynamics of TIs are applicable in the low-energy limit, optical experiments should be performed at terahertz (THz) or sub-THz frequencies [9][10][11][12] . Qualitatively, equation (1) applied to the magnetic and electric fields of the primary THz radiation results in a perpendicular polarized secondary THz radiation. The sum of the primary and secondary radiation can be viewed as the rotation of the polarization plane, that is, as the Faraday effect. We would like to note that the quantum Faraday effect and the TME are basically different manifestations of the same axion physics 13 .
In real samples, the TME may be screened by nontopological contributions [13][14][15] . In fact, quantized Faraday rotation has been detected in 2D electron gas 16 and graphene 17 in the quantum Hall effect (QHE) regime. However, in both experiments the Faraday rotation takes a value y F ¼ À 4a/(1 þ n sub ), that is, it depends on the refractive index of the substrate n sub and hence is not fundamental.
Here we report on the observation of the universal Faraday rotation angle equal to the fine structure constant a. Strained HgTe layers grown on CdTe, that are investigated in the present work, are shown to be a 3D TI 18 with surface-dominated charge transport 19 that was observed in THz experiments as well 9,11,20 . To eliminate the material details, we perform measurements under antireflection conditions, such that the transmission through the CdTe substrate is approaching 100% (ref. 21).

Results
THz Faraday effect in 2D systems. The observed Faraday rotation angle y F j j¼ a (for N ¼ 0) comes from two spatially separated topological surfaces in a 3D TI. This corresponds to the half-quantized Hall conductivity e 2 /(2h) per surface or, equivalently, to the TME occurring at each surface separately. Therefore, the observed Faraday effect 2(N þ 1/2)a is intimately related to the TME, which distinguishes qualitatively our 3D TI from 2D or quasi-2D materials. There is also a quantitative difference. Even without the substrate (n sub ¼ 1), the Faraday rotation in graphene would be quantized as 4(N þ 1/2)a, including the spin and valley degeneracies. The minimum Faraday rotation angle is then y F j j¼ 2a (for N ¼ 0) 17 . For a conventional 2D electron gas, such as hosted in GaAs/AlGaAs heterostructures, the quantization of the Faraday angle is expected to be 2Na, where the factor of 2 comes from the equal contributions of the up-and down-spin subsystems, which independently exhibit the integer QHE. This is because in GaAs/AlGaAs heterostructures, the Zeeman splitting for magnetic fields below 10T is negligible compared to the THz photon energy. Therefore, the minimum Faraday rotation angle would be also y F j j¼2a (for N ¼ 1) 16 , which is twice larger than our result.
Sample details. The strained HgTe film is a 58 nm thick HgTe layer embedded between two Cd 0.7 Hg 0.3 Te layers (Fig. 1a). The Cd 0.7 Hg 0.3 Te layers have a thickness of 51 nm (lower layer) and 11 nm (top/cap layer), respectively. The purpose of these layers is to provide the identical crystalline interface for top and bottom surface of the HgTe films as well as to protect the HgTe from oxidization and adsorption. This leads to an increase in carrier mobility with a simultaneous decrease in carrier density compared to uncaped samples 18 . The transport characterization on a standard Hall bar sample shows a carrier density at 0 V gate of 1.7 Â 10 11 cm À 2 and a carrier mobility of 2.2 Â 10 5 cm 2 V À 1 s À 1 . The optical measurements are carried out on a sample fitted with a 110 nm thick multilayer insulator of SiO 2 /Si 3 N 4 and a 4 nm thick Ru film. The Ru film (oxidized in the air) is used as a semitransparent top-gate electrode 22 . The gate leads to B15% suppression of the transmission signal, which can be taken into account as field-independent contribution to the conductivity. The properties of the gate material have been investigated in a separate experiment.
THz spectroscopy. The transmittance experiments at THz frequencies (0.1 THzovo1 THz) are carried out in a Mach-Zehnder interferometer arrangement 23,24 , allowing measurement of the amplitude and the phase shift of the electromagnetic radiation in a geometry with controlled polarization (Fig. 1a). The monochromatic THz radiation is provided by a backward-wave oscillator. The THz power on the sample is in between 10 and 100 mW with the focal spot of 0.2 cm 2 . Using wire grid polarizers, the complex transmission coefficient t ¼ t j je if is obtained both in parallel t p (Fig. 1b) and cross t c (Fig. 1c) polarization geometries, providing full information about the transmitted light. External magnetic fields B 7 T are applied using a split-coil superconducting magnet. The experiments are carried out in Faraday geometry, that is, with B applied parallel to the propagation direction of the THz radiation. The ac conductivity tensorŝðoÞ at THz angular frequency o ¼ 2pn is obtained from the experimental data by inverting the Berreman equations 25 for the complex transmission coefficient through a thin conducting film on an insulating substrate. The explicit expressions used in the calculations are given in Methods section.
In general case, the light propagating along the z direction can be characterized by the orthogonal x and y components of the electric and magnetic fields, which can be written in the form of a 4D vector V. The interconnection between vectors V 1 and V 2 , corresponding to different points in space separated by a distance ', is given by V 1 ¼Mð'ÞV 2 . HereMð'Þ is a 4 Â 4 transfer matrix. For an insulating substrate of thickness ' and dielectric constant e, this is the identity matrixM CdTe ð'Þ¼I provided 2' ffiffi e p n=c is an integer. We find in a separate experiment on a bare CdTe substrate that this condition is fulfilled for vE0. 35 THz, and all the measurements presented here are performed at this frequency to minimize the contribution to the Faraday signal from the substrate. The corresponding photon energy of 1.4 meV is much smaller than the energy gap in strained HgTe (above 10 meV) 18 , and equations (1) are a good approximation.
For normal incidence, the fields across the conducting interface are connected by the Maxwell equation rÂH¼ŝE. Here the e À iot time dependence is assumed for all fields. As the wavelength of 856 mm for v ¼ 0.35 THz is much larger than the HgTe layer thickness, we use the limit of thin film, and the corresponding transfer matrixM HgTeŝ ð Þ is determined by the diagonal (s xx ) and Hall (s xy ) components of the conductivity tensorŝ. Within the Drude-like model, these components for one type of charge carriers can be written in the form 13,26 Here O c is the cyclotron resonance (CR) frequency, s 0 is the dc conductivity, and t is the scattering time. For classical conductors, the CR frequency is written where m e is the effective electron cyclotron mass. The total transfer matrixM¼M CdTeMHgTe relates vectors V on both sides of the sample and hence contains full information about the transmission and reflection coefficients. Thus, when M CdTe is the identity matrix, the influence of the substrate is minimized, and the THz response is dominated by the ac transport properties of the HgTe layer, in accord with equations (2) and (3). The calculation of the complex transmission coefficients t p and t c based on the transfer matrix formalism as well as the exact form of the transfer matrices are presented in Methods section.
Magnetic field dependence of the THz transmission is dominated by a sharp CR of surface electrons (e) O ce at B e ¼0:4 T (Fig. 1b,c). Below we demonstrate their Dirac-like character and that they are responsible for the universal Faraday rotation. Remarkably, the observation of the CR both in t p and t c indicates a high purity of our HgTe layer. The scattering time is significantly longer than the inverse THz frequency ot ) 1, and according to equations (2) and (3) the ac conductivity reveals a resonance-like behaviour s xx , s xy p1/(O 2 ce À o 2 ). Further features are broad resonances at B s1 ¼3:7 T and at B s2 ¼2:2 T indicated in Fig. 1 as s1 and s2, respectively. The phase of the corresponding THz transmission coefficient f c in the vicinity of these resonances has the opposite sign with respect to that of the e-CR. Remarkably, the s1 and s2 resonances disappear with applying positive gate voltage (Figs 1 and 2b). We associate them with either interband Landau-level transitions or thermally activated states as discussed below.
Band structure analysis. To understand the origin of the experimentally observed resonances, we analyse the band structure of tensile strained Cd 0.7 Hg 0.3 Te/HgTe layer as shown in Fig. 2a. It is obtained similar to ref. 19 within the tight binding approximation of the 6 Â 6-Kane Hamiltonian 27,28 . Due to reduced point symmetry at the boundary between the Cd 0.7 Hg 0.3 Te and HgTe layers, an additional interface potential is allowed in the Hamiltonian 29 . This potential is used to shift the Dirac point closer to the valence band edge, so that the tight binding results are in good agreement with recent angle-resolved photoemission spectroscopy (ARPES) experiments 18,30 and ab initio calculations 31 on HgTe. Final position of the Dirac point is at energy around À 40 meV (that is, it is burried in the the valence band) in agreement with 18,30,31 . The Dirac-like surface states are located in the band gap between the light-hole (LH, conduction) and heavy-hole (HH, valence) subbands (see the red line in Fig. 2a), and they are not perfectly linear due to the hybridization with the heavy-hole band. The dashed red line in Fig. 2a shows the linear dispersion of the 2D surface state without hybridization with the heavy-hole band. This is a hypothetical curve, since in the realistic material the hybridization with the heavy-hole band is non-zero. The camel back of the heavy-hole band originates from coupling of this band to the electron-like valence band and is therefore a hallmark of the inverted band structure of HgTe. In accordance with previous transport data, the chemical potential crosses the topological surface states for a large range of gate voltages 19 . The total electron density in Fig. 2a is n % 2Â10 11 cm À 2 , representing the experimental situation at U G ¼1:9 V. For simplicity, we assume here the same density at the top and bottom surfaces. Using the general formula for a quasiclassical CR 32 is the energy dispersion, B is the magnetic field and A is the area enclosed by the wave vector k, we calculate for the topological surface state O ce =B % 35 cm À 1 T À 1 .
Quantized THz Hall effect. Experimentally, simultaneous fit of the real and imaginary parts of t p and t c allows the extraction of all transport characteristics, that is, conductivity, charge carrier density, scattering time and CR frequency 21 . The inset of Fig. 2b shows experimentally determined electron CR as a function of gate voltage, which perfectly agrees with the theoretical value for the topological Dirac-like surface states. Since only surface states are observed in transport experiments on the similar structures 19 , a possible explanation of the appearance of additional resonances is interband Landau-level transitions between the HH bulk bands and topological surface states. Such transitions are generally allowed as can be shown using the Kubo formula. Another possibility would be thermally activated transport between the camel back of the HH bulk band and the surface states. This is generally possible since the THz field may well induce heating of the carriers, resulting in a higher effective temperature compared to that of the lattice.
From the obtained scattering time and the CR positions in the magnetooptical spectra of Fig. 1b,c, one can calculate the mobility m ¼ tO c /B. The surface states demonstrate high mobility m e ¼1:8Â10 5 cm 2 V À 1 s À 1 , which agrees with the dc transport data. Since the e-CR and s1,s2 resonances occur at different magnetic fields, their contributions to the ac transport can be clearly separated, as presented in Fig. 2b. The striking feature of this plot is that the ac conductivity of the surface states dominates at large gate voltages. In what follows, we concentrate therefore on U G 41:0 V, while remaining weak contribution from the interband Landau-level transitions/thermally activated transitions are subtracted as explained in Supplementary Fig. 1. Figure 3a demonstrates the real part of the Hall conductivity s xy . The overall behaviour is provided by the high-field tail of the classical Drude model, that is, equation (3), resulting in a rapid suppression of s xy with growing magnetic field. In addition to the classical behaviour, regular oscillations in @s xy /@B can be recognized, which are linear in inverse magnetic field (Fig. 3b). The slope of the linear behaviour changes with gate voltage, reflecting gate dependence of the electron density per surface. These QHE oscillations are extrapolated to 1/2, indicating Dirac character of the surface electrons 33 . The oscillations of @s xy /@B in Fig. 3a are superimposed by the 1/B tail of the classical electron CR and therefore are not well pronounced. The visibility can be significantly improved by inserting the sample in a Fabry-Pérot resonator, as we have previously demonstrated for a similar structure 11 . In magnetic fields above 5 T, the Hall conductivity clearly shows a plateau close to s xy ¼ e 2 /h, corresponding to a value (1/2)e 2 /h per surface (Fig. 3a). Another plateau close to (3/2)e 2 /h per surface is also recognizable at a magnetic field of 3 T. It is superimposed on the CR tail and therefore tilted. The steps in s xy loose their regularity in lower-magnetic fields, as can be qualitatively explained by a finite THz frequency o in magnetooptical experiments. As mentioned above, the overall behaviour of s xy (o) is provided by the classical curve of equation (3), and the real part of s xy can be approximated as s xy (o)Es 0 O ce /[(O 2 ce À o 2 )t], which in the limit O ce ) o reduces to the expression s xy ¼ ne/B, being a multiple of e 2 /h. In low magnetic fields, the CR frequency O ce becomes comparable to the THz frequency o, destroying the regularities in s xy (o).
Since in strained HgTe, the Fermi level lies in the bulk band gap (see Fig. 2a and ref. 19), we attribute the observed THz QHE to the formation of the 2D Landau levels at the top and bottom surfaces of the HgTe layer (Fig. 1a). This interpretation is further substantiated by our theoretical analysis of the ac quantum Hall conductivity s xy (o) calculated from the Kubo formula for both top and bottom surface states within the Dirac model 13,34 .
Our two-surface Dirac model describes well the surface carrier CR (the inset of Fig. 3a). The lengths of the theoretical Hall plateaus in the high magnetic field region (Fig. 3a) correlate correctly with the positions of the extrema in the derivative @s xy /@B. However, the model predicts much sharper transitions between the QHE plateaus, as observed in the experiment. One of two possible explanations is the heating of the surface carriers by the THz field, resulting in a higher effective temperature compared to that of the lattice. Such a heating can occur due to inefficient energy relaxation in the electronic system through the emission of LO phonons at low temperatures 35 . The best fit of our experimental data is obtained with T ¼ 25 K (Fig. 3a). Another explanation is based on spatial fluctuations of the surface carrier densities, which are likely to occur in our samples due to their large lateral sizes compared to the typical Hall bars used in the dc measurements. The experimental data of Fig. 3a can alternatively be well fitted assuming cold carriers (T ¼ 1.8 K) with density fluctuations within 10% relative to their nominal values. As the fits are nearly indistinguishable, we cannot quantitatively determine the contributions of both mechanisms leading to the smearing of the THz QHE plateaus.
The stronger the field, the closer the Hall conductivity to the quantized values expected for a two-surface Dirac system where N a,b are the integer numbers of the highest occupied Landau levels at the top and bottom surfaces, with n a En b being the corresponding carrier densities (F 0 ¼ h/ e j j is the magnetic flux quantum). Upon approaching the CR, the Hall conductivity deviates from the quantized values in equation (4) due to the predominance of the intraband transitions between the Landau levels. From the fitting procedure, we extract the nominal carrier densities n a ¼0:92Â10 11 cm À 2 and n b ¼1:07Â10 11 cm À 2 . The total surface carrier density n a þ n b agrees well with that obtained from Drude-like fits of magnetooptical spectra. Another fitting parameter is the classical (Drude) surface conductivity s a ¼ s b E50e 2 /h. Its large value indicates high-surface carrier mobility, insuring that the condition for the quantum Hall Þ 1=2 is the characteristic Landau-level spacing for a Dirac system, t a,b are the scattering times of the top and bottom carriers, R 0 ¼ h/(2e 2 ) is the resistance quantum, and v F is the Fermi velocity of the Dirac surface states.
Quantized THz Faraday effect. Having established that the THz response of the topological surface states in high magnetic fields B45 T is determined by the conductivity quantum G 0 ¼ e 2 /h (N a ¼ N b ¼ 0), we turn to the central result of this work, the THz Faraday effect. Owing to the TME of equation (1), an oscillating electric E x e À iot (magnetic H y e À iot ) field of the linearly polarized THz radiation induces in a 3D TI an oscillating magnetic aH x e À iot (electric aE y e À iot ) field. The generated, in such a way, secondary THz radiation is polarized perpendicular to the primary polarization and its amplitude is a times smaller. This can be viewed as a rotation of the initial polarization by an angle y F j j¼ arctan a % 7:3Â10 À 3 rad. Indeed, Fig. 4a clearly demonstrates that the Faraday angle in high magnetic fields is close to this fundamental value.
We rigorously characterize the THz Faraday effect, and the Faraday ellipticity Z F is shown in Fig. 4b. It is relatively small Z F j jo y F j j in high magnetic fields, but does not reach zero. This observation indicates that while the TME dominates, the interaction of TIs with THz radiation is not a completely dissipationless process in our samples. Remarkably, the universal value of the Faraday angle remains robust against the gate voltage. This is demonstrated in Fig. 4c, where y F j jEa for 0:7 VoU G o1:9 V.

Discussion
The observed terahertz Faraday rotation equal to the fine structure constant a ¼ e 2 /2E 0 hc ¼ e 2 m 0 c/2h is a direct consequence of the TME, confirming axion electrodynamics of 3D topological insulators. We use monochromatic terahertz spectroscopy, providing complete amplitude and phase reconstruction, which can be applied to investigate topological phenomena in various systems, including graphene, 2D electron gas, layered superconductors and recently experimentally discovered Weyl semimetals 36 . Picoradian angle resolution can be achieved using a balanced detection scheme 37 , and the universal Faraday rotation in combination with the magnetic flux quantum F 0 ¼ h/ e j j and the conductivity quantum G 0 ¼ e 2 /h are suggested 14 to use for a metrological definition of the three basic physical constants, e, h and c (given that the vacuum permeability is equal exactly m 0 ¼4pÂ10 À 7 N=A 2 ).

Methods
Experimental technique. Magnetooptical experiments in the THz frequency range (100 GHzono1,000 GHz) have been carried out in a Mach-Zehnder interferometer arrangement 23 , which allows measurements of the amplitude and the phase shift in a geometry with controlled polarization of radiation.
Theoretical analysis of magnetooptical spectra. To analyse the experimental transmission spectra, we follow the formalism described by Berreman 24,25 . The THz light propagating along the z direction can be characterized by the (orthogonal) components of electric (E x , E y ) and magnetic (H x , H y ) fields, which may be combined in form of a four-component vector V ¼ (E x , E y , H x , H y ). The propagation of light between two points in space separated by a distance ' and characterized by vectors V 1 and V 2 can be described via 4 Â 4 transfer matrixMð'Þ as V 2 ¼Mð'ÞV 1 . To provide a simple example, for an isotropic dielectric substrate  Fig. 1).
the transfer matrix is simplified to: Here Z ¼ ffiffiffiffiffiffiffi m=e p and k¼ ffiffiffiffiffi me p o=c is the wave vector. In the following, we assume m ¼ 1. The Berreman procedure is in general not limited to the case of normal incidence. However, in such geometry the choice of tangential field components simplifies the treatment. Electric and magnetic fields across the interfaces are connected by the Maxwell equation rÂH¼ŝE. Hereŝ is the complex conductivity tensor of the material and the time dependence in form e À iot is assumed. The full transfer matrixM¼M CdTeMHgTe describes the transmission and reflection coefficients, which can be calculated using another basis. In this basis, the vector V consists of (i) the amplitude of the linearly polarized wave (E x ) propagating in the positive direction, (ii) the amplitude of the wave with the same polarization propagating in the negative direction, and (iii) and (iv) of two waves with perpendicular polarization (E y ). The propagation matrix in the new basis isM 0 ¼V À 1MV , witĥ The present experiment is described by a linearly polarized incident wave and by two components of the transmitted (t) and reflected (r) waves, respectively. The equation connecting all waves is given by: Here the t p and t c are the complex transmittance amplitudes within parallel and crossed polarizers geometry, r p and r c are respective reflectivity coefficients. The Faraday rotation y and Faraday ellipticity Z are obtained from the transmission amplitudes t p , t c j j and phase shifts f p , f c as To interpret the experimental data, we use the ac conductivity tensorŝðoÞ obtained in the classical (Drude) limit from the Kubo conductivity of topological surface states 7,15 . We note that for a 2D conducting film on an isotropic dielectric substrate, the complex transmission can be obtained analytically. For a substrate with permittivity e, the final equations for the spectra in parallel (t p ) and crossed (t c ) polarizers are given by: where a xy ¼s xy Z 0 cosðk'Þ À iZ sinðk'Þ ð Þ : Here ' is equal to the substrate thickness, Z 0 E377 Ohm is the impedance of free space, is the relative impedance of the substrate and k ¼ ffiffi e p o=c is the wave vector of the radiation in the substrate. The components of the conductivity tensor, s xx (o), and s xy (o), are given by equations (2) and (3). As can be clearly seen, equations (9) and (10) can be inverted analytically to get the solution for the complex conductivity thus avoiding the numerical procedure.
As has been discussed previously 6 , in case of a film on a substrate the universal value of the Faraday rotation angle should be modified by the refractive index of the substrate. We note that an infinite substrate has been assumed in these calculations. In the present case of finite substrate, exact transmission matrix formalism has been utilized thus automatically taking into account the influence of the substrate and of the Ru-gate. Moreover, monochromatic radiation has been used in the experiments and the frequency of this radiation has been selected close to a maximum of the Fabry-Pérot resonances in the substrate (sinðk'Þ ! 0). In this case, equations (11) and (12) can be simplified to a xx ¼2 þ s xx Z 0 ð13Þ a xy ¼s xy Z 0 ð14Þ and the influence of the substrate is minimized.
Data availability. The data that support the findings of this study are available from the corresponding authors upon reasonable request.