Reversed electromagnetic Vavilov-\v{C}erenkov radiation in naturally existing magnetoelectric media

TL;DR

This paper theoretically demonstrates the possibility of reversed Vavilov-Cerenkov radiation in naturally occurring magnetoelectric media, especially topological insulators, by analyzing electromagnetic responses at interfaces and particle motion.

Contribution

It introduces a Green's function approach for magnetoelectric interfaces and predicts reversed Cerenkov radiation in natural materials like topological insulators.

Findings

Reversed Vavilov-Cerenkov radiation can occur in magnetoelectric media.

Estimated frequencies of reversed radiation are comparable to experimental values in metamaterials.

Green's function method can be generalized to various geometries.

Abstract

We consider two semi-infinite magnetoelecteric media separated by a planar interface whose electromagnetic response is described by axion-electrodynamics. The time dependent Green's function characterizing this geometry is obtained by a method which can be directly generalized to cylindrical and spherical configurations of two magnetoelectrics separated by an interface. We establish the far-field approximation of the Green's function and apply these results to the case of a charged particle moving from one medium to the other at a high constant velocity perpendicular to the interface. From the resulting angular distribution of the radiated energy per unit frequency we provide theoretical evidence for the emergence of reversed Vavilov-\v{C}erenkov radiation in naturally existing magnetoelectric media. In the case where one of the magnetolectrics is a 3D topological insulator, ${\rm TlBiSe}_2$ for example, located in front of a regular insulator, we estimate that an average forward Vavilov-\v{C}erenkov radiation with frequency $\sim 2.5 \,\, {\rm eV}$ ($\sim 500\,\, {\rm nm}$) will produce a highly suppressed reversed Vavilov-\v{C}erenkov radiation which can be characterized by an effective frequency in the range of $\sim (4\times 10^{-3}-0.5) \,\, {\rm meV}$. However, this value compares favorably with recent measurements in left-handed metamaterials yielding reversed Vavilov-\v{C}erenkov radiation with frequencies of the order of $(1.2-3.9)\times 10^{-2}\,\, {\rm meV}$.

Tables

Table 1. Table 1: This table shows the orders of magnitude for the power radiated per unit frequency in the backward direction (reversed VC radiation) when a particle with v = 0.8 𝑣 0.8 v=0.8 and charge q = α 𝑞 𝛼 q=\sqrt{\alpha} propagates across the interface of a normal insulator and the topological insulator TlBiSe 2 subscript TlBiSe 2 {\rm TlBiSe}_{2} characterized by n = 2 𝑛 2 n=2 , θ ~ = 11 ​ α ~ 𝜃 11 𝛼 \tilde{\theta}=11\alpha .

$\omega \left(\mathrm{eV}\right)$	$d^2{E}_{\mathrm{RVCR}}/dtd\omega$ $\left(\mu W/\mathrm{eV}\right)$
$2$	$3.5\times {10}^{-3}$
$4$	$7.0\times {10}^{-3}$
$6$	$1.0\times {10}^{-2}$
$8$	$1.4\times {10}^{-2}$

Equations

$${\alpha }_{ij}={\left[\frac{\partial M_j}{\partial E_i}\right]}_{\mathbf{B}=0}={\left[\frac{\partial P_i}{\partial B_j}\right]}_{\mathbf{E}=0}\equiv {\tilde{\alpha }}_{ij}+\frac{e^2}{\hslash c}\frac{\vartheta }{4{\pi }^2}{\delta }_{ij},$$

$$\frac{S_{\vartheta }}{\hslash }=\frac{\vartheta }{32{\pi }^2}\int d^{4}x{ϵ}^{\alpha \beta \mu \nu }\frac{1}{e^2}{F}_{\alpha \beta }{F}_{\mu \nu }=\vartheta C_2,$$

$${\sigma }_H^{\Sigma }=\frac{e^2}{h}\left(\left|\frac{1}{2}\right|+m\right),$$

$$\vartheta (z)=H(z){\vartheta }_2+H(-z){\vartheta }_1,$$

$$S[\Phi ,\mathbf{A}]$$

$$\mathbf{E}=-\frac{\partial \mathbf{A}}{\partial t}-\nabla \Phi ,\mathbf{B}=\nabla \times \mathbf{A},$$

$$\nabla \cdot \mathbf{B}=0,\nabla \times \mathbf{E}=-\frac{\partial \mathbf{B}}{\partial t}.$$

$$ϵ\nabla \cdot \mathbf{E}$$

$$\nabla \times \mathbf{B}-ϵ\frac{\partial \mathbf{E}}{\partial t}$$

$$\tilde{\theta }=\alpha ({\vartheta }_2-{\vartheta }_1)/\pi .$$

$$\tilde{\theta }=\alpha (2m+1),$$

$${\varrho }_{\vartheta }=\frac{1}{4\pi }\tilde{\theta }\delta (z)\mathbf{B}\cdot \hat{\mathbf{u}},{\mathbf{J}}_{\vartheta }=\frac{1}{4\pi }\tilde{\theta }\delta (z)\mathbf{E}\times \hat{\mathbf{u}},$$

$$\nabla \cdot {\mathbf{J}}_{\vartheta }+\frac{\partial {\varrho }_{\vartheta }}{\partial t}=0$$

$$\left[{\partial }_i\delta (z)\right]{\epsilon }_{ijk}{E}_j{u}_k={\delta }^{\prime }(z){\delta }_i^3{\epsilon }_{ijk}{E}_j{u}_k=0,$$

$$\nabla \cdot \mathbf{D}=4\pi \varrho$$

$$\nabla \cdot \mathbf{B}=0$$

$$\mathbf{D}=ϵ\mathbf{E}-\frac{\alpha }{\pi }\vartheta (z)\mathbf{B},\mathbf{H}=\mathbf{B}+\frac{\alpha }{\pi }\vartheta (z)\mathbf{E}.$$

$$ϵ{\left[{\mathbf{E}}_z\right]}_{z=0^{-}}^{z=0^{+}}=\tilde{\theta }{\mathbf{B}}_z{|}_{z=0},{\left[\hat{\mathbf{u}}\times \mathbf{B}\right]}_{z=0^{-}}^{z=0^{+}}=-\tilde{\theta }(\hat{\mathbf{u}}\times \mathbf{E}){|}_{z=0},$$

$${\left[{\mathbf{B}}_z\right]}_{z=0^{-}}^{z=0^{+}}=0,{\left[\hat{\mathbf{u}}\times \mathbf{E})\right]}_{z=0^{-}}^{z=0^{+}}=0,$$

$$\left[{\left[{\square }^2\right]}_{\nu }^{\mu }-\tilde{\theta }\delta (z){\epsilon }_{\nu }^{3\mu \alpha }{\partial }_{\alpha }\right]A^{\nu }=4\pi J^{\mu },$$

$$ϵ\frac{\partial \Phi }{\partial t}+\nabla \cdot \mathbf{A}=0,$$

$${\left[{\square }^2\right]}_{\nu }^{\mu }=(ϵ{\square }^2,{\square }^2{\delta }_j^i),{\square }^2=ϵ{\partial }_t^2-{\nabla }^2.$$

$$A^{\mu }{|}_{z=0^{-}}^{z=0^{+}}=0$$

$${\partial }_z{A}^i{|}_{z=0^{-}}^{z=0^{+}}$$

$$\left[{\left[{\square }^2\right]}_{\nu }^{\mu }-\tilde{\theta }\delta (z){\epsilon }_{\nu }^{3\mu \alpha }{\partial }_{\alpha }\right]G_{\sigma }^{\nu }(x,x^{\prime })=4\pi {\eta }_{\sigma }^{\mu }{\delta }^4(x-x^{\prime }),$$

$$A^{\mu }(x)=\int d^4{x}^{\prime }{G}_{\nu }^{\mu }(x,x^{\prime })J^{\nu }(x^{\prime }),$$

$$G_{\nu }^{\mu }(x,x^{\prime })$$

$$\left[{\mathcal{O}}^{\mu }{}_{\nu }+i\tilde{\theta }\delta (z){\epsilon }_{\nu }^{3\mu \alpha }{k}_{\alpha }\right]g_{\sigma }^{\nu }(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )={\eta }_{\sigma }^{\mu }\delta (z-z^{\prime }),$$

$${\mathcal{O}}^{\mu }{}_{\rho }{\mathcal{H}}_{\nu }^{\rho }(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )={\eta }_{\nu }^{\mu }\delta (z-z^{\prime })$$

$$ϵ\left({\mathbf{k}}_{\perp }^2-{\omega }^2ϵ-{\partial }_z^2\right){\mathcal{H}}_0^0(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )$$

$$ϵ\left({\mathbf{k}}_{\perp }^2-{\omega }^2ϵ-{\partial }_z^2\right){\mathcal{H}}_i^0(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )$$

$$({\mathbf{k}}_{\perp }^2-{\omega }^2ϵ-{\partial }_z^2){\mathcal{H}}_0^i(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )$$

$$({\mathbf{k}}_{\perp }^2-{\omega }^2ϵ-{\partial }_z^2){\mathcal{H}}_j^i(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )$$

$${\mathcal{H}}_i^0(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )=0={\mathcal{H}}_0^i(z,z^{\prime };{\mathbf{k}}_{\perp },\omega )$$

Full text

I Introduction

Since its experimental discovery in $1934$ cherenkov ; Vavilov , Vavilov-Čerenkov (VC) radiation has played a special role in the study of high-energy particle physics, high-power microwave sources and nuclear and cosmic-ray physics Jelley ; Jelley1 , both theoretically and phenomenologically. VC radiation occurs when charged particles travel across a dielectric medium and propagate through the material with velocity $v$ higher than $c/\sqrt{\epsilon\mu}$. Here $c$ is the speed of light in vacuum, $\epsilon$ is the permittivity of the medium, $\mu$ is its permeability, which we take equal to one, and $n=\sqrt{\epsilon}$ is the refraction index of the material. Throughout this paper, we use Gaussian units. The first theoretical description of such radiation in the framework of Maxwell’s theory, developed by Frank and Tamm in Ref. FT , revealed its unique polarization and directional properties. In particular, VC radiation is produced in a forward cone defined by the angle $\theta_{C}=\cos^{-1}[c/vn]<\pi/2$ with respect to the direction of the incident charge. Since the emergence of accelerators in nuclear and high-energy physics, VC radiation has been widely used to design an impressive variety of detectors, such as e.g. the ring imaging Čerenkov detectors Ypsilantis , which can identify charged particles by providing a straightforward effective tool to test its physical properties like velocity, energy, direction of motion and charge PADG . As remarkable cases, the antiproton Chamberlain and the J particle Aubert were discovered using Čerenkov detectors.

In recent years, the study of the reversed VC radiation, an exotic electromagnetic phenomenon, has long been of considerable interest Pendry ; Pendry2 ; LU ; Luo ; RRCR ; Sheng ; Chen ; Duan ; Tao . The reversed VC radiation occurs when the photons are emitted in the backward direction with respect to the velocity of the propagating charged particle. This remarkable analytic prediction was introduced by Veselago in Ref.Veselago , invoking materials having simultaneously a negative permittivity and permeability, dubbed as left-handed media (LHM), as opposed to normal dielectrics called right-handed media (RHM). This distinction is better characterized in terms of an electromagnetic wave with electric field ${\mathbf{\bar{E}}}$, magnetic field ${\mathbf{\bar{B}}}$, and wave vector ${\mathbf{k}}$ propagating in a linear medium where the conditions ${\mathbf{k}}\cdot{\mathbf{\bar{E}}}=0$ and ${\mathbf{k}}\cdot{\mathbf{\bar{B}}}=0$ hold. A medium is said to be left-handed (right-handed) according to ${\mathbf{k}}$ being antiparallel (parallel) to ${\mathbf{\bar{E}}}\times{\mathbf{\bar{B}}}$, i. e., when the group velocity and the phase velocity of the wave are antiparallel (parallel). The realization of a negative refraction index by creating an interface between a LHM and an ordinary media was also considered in Ref. Veselago and experimentally proven later in Ref. SHELBY . Since LHM are not readily available in nature, they have been artificially constructed in the laboratory as metamaterials, e.g., by combining metallic thin wires Pendry with split-ring resonators Pendry2 in a periodic cell structure. Also plasmonic thin films LU and photonic crystals Luo constitute an alternative.

In this paper, we show that reversed VC radiation occurs in magnetoelectric media ODELL , which are naturally existing RHM like antiferromagnets Hehl , topological insulators (TIs) TI1 ; TI11 ; QILIBRO , and Weyl semimetals ARM , for example. It is interesting to observe that Reversed Čerenkov sound in topological insulators has been already discussed in Ref. SMIRNOV . Linear magnetoelectric media are characterized by the magnetoelectric polarizability tensor $\alpha_{ij}$, which is defined as the response of the magnetization under a change of an electric field or, equivalently, as the response of the polarization under a change of a magnetic field, i.e.,

[TABLE]

where $\vartheta$ is a dimensionless parameter. The equality follows because $M_{j}=\partial{\mathcal{E}}/\partial B_{j}$ and $P_{i}=\partial{\mathcal{E}}/\partial E_{i}$, with ${\mathcal{E}}$ being the electromagnetic enthalpy density VDBPRL102_2009 . Here ${\tilde{\alpha}}_{ij}$ denotes the traceless part of $\alpha_{ij}$ . Leaving aside the remarkable microscopic properties of magnetoelectric media, we will consider only the effective-medium theory describing their electromagnetic response arising from the term proportional to $\delta_{ij}$ in Eq. (1), i. e., assuming an isotropic material. Such effective description is provided by the standard Maxwell Lagrangian ${\cal L}_{\rm ED}$ plus the additional term ${\cal L}_{\vartheta}=-\frac{e^{2}}{\hbar c}\frac{\vartheta}{4\pi^{2}}\,\mathbf{E}\cdot\mathbf{B}$. Here $\mathbf{E}$ and $\mathbf{B}$ are the electromagnetic fields, $\alpha$ is the fine-structure constant and $\vartheta$ is a field known as the (scalar) magnetoelectric polarizability (MEP) in condensed matter physics or the axion field in particle physics axion . In this way, the addition of ${\cal L}_{\vartheta}$ to the Maxwell Lagrangian is usually referred to as the Lagrangian density for axion electrodynamics. However, we consider $\vartheta$ as an additional parameter characterizing the material, in the same footing as its permittivity $\epsilon$ and permeability $\mu$, which lead us to restrict the designation of axion electrodynamics to $\vartheta$ electrodynamics ($\vartheta$-ED) emphasizing that $\vartheta$ is not a dynamical field. The nature of the MEP depends on the type of magnetoelectric material under consideration and it is deeply related to the magnetic symmetries of the substance DZYA ; Rivera and/or to the properties of its band structure TI11 . It can be calculated from a Kubo-type response formula, once a microscopic model Hamiltonian for the material is adopted. The permittivity $\epsilon$ is usually designed by the Drude-Lorentz type of single resonance oscillator model Mohideen . As as first step in dealing with radiation, we consider $\vartheta$ and $\epsilon$ frequency independent.

The main feature characterizing magnetoelectric materials, which is responsible for most of their unusual effects, is the magnetoelectric (ME) effect arising from the additional contribution ${\cal L}_{\vartheta}$ QI_SCIENCE . This coupling produces effective field-dependent charges and current densities, which allow the generation of an electric (magnetic) polarization due to the presence of a magnetic (electric) field. Let us emphasize that we are dealing with standard electrodynamics of a RHM ($\epsilon>1,\,\mu>1,\,n>1$) supplemented only by additional field-dependent sources.

Among the choices of accessible magnetoelectric media, we consider materials characterized by a piecewise constant value of the MEP $\vartheta$ and focus on (TIs). Three-dimensional (3D) strong TIs are a fascinating class of RHM that can host a conducting helical surface state with an odd number of fermions at the low-energy limit, each having the dispersion relation of a nondegenerate Dirac cone with a crossing point at or close to the Fermi level $E_{\rm F}$. These gapless crossing points are called Dirac points. Nevertheless, TIs behave as magnetoelectric insulators in the bulk with a finite energy gap. The surface state is further topologically protected by time-reversal symmetry and/or inversion symmetry, coupled with spin-momentum locking properties. The latter means that the spin orientation of the electrons on the surface Dirac cone is always locked perpendicularly to their momentum. A distinguishing feature of 3D strong TIs among magnetoelectrics, is that the dimensionless parameter $\vartheta$, dubbed as the orbital magnetoelectric polarizability in this case, is of topological nature and arises from the bulk band structure. It is given by a non-Abelian Berry flux over the Brillouin zone and results in an integer multiple of $\pi$ TI11 ; VDBPRL102_2009 ; QZh .

Some general comments regarding the properties of the additional term ${\cal L}_{\vartheta}$ in the case of TIs are now in order. Let us recall that the electric and magnetic fields have dimensions of charge divided by distance squared in Gaussian units (for the moment we retain $\hbar$ and $c$). In this way, the contribution of the Lagrangian density ${\cal L}_{\vartheta}$ to the action is $S_{\vartheta}=(c\hbar^{2}/e^{4})\int dt\,d^{3}x\,{\cal L}_{\vartheta}$. Rewriting ${\mathbf{E}}\cdot{\mathbf{B}}$ in terms of the field strength tensor $F_{\alpha\beta}=\partial_{\alpha}A_{\beta}-\partial_{\beta}A_{\alpha}$, where $A_{\mu}$ is the electromagnetic potential, and considering a closed spacetime with no boundaries, we get

[TABLE]

where $C_{2}$ is an integer. This is because in such spaces the dimensionless integral in Eq. (2) is equal to $32\pi^{2}\,C_{2}$, where $C_{2}$ is the second Chern number of the manifold FUJIKAWA . Under changes of $\vartheta$, the quantity $\exp(-iS_{\vartheta}/\hbar)$ must remain invariant, which means that two values of $\vartheta$ differing by an integer multiple of $2\pi$ are equivalent. Further imposing time-reversal (TR) symmetry yields new constraints on the values of $\vartheta$. Since $\mathbf{E}\cdot\mathbf{B}$ is odd under TR, one could think that the only allowed value would be $\vartheta=0$ (modulo $2\pi)$. Nevertheless, the condition $\exp(-iS_{\vartheta}/\hbar)=\exp(+iS_{\vartheta}/\hbar)$ yields the possibility of having $\vartheta=\pi$. In this way, we obtain two families of magnetoelectric materials described by the choices $\vartheta_{1}=0$ (normal insulators) and $\vartheta_{2}=\pi$ (topological insulators). Both values of $\vartheta$ are defined modulo $2\pi$.

Another important property of $S_{\vartheta}$ is that the integrand in Eq. (2) is a total derivative, as can be readily seen by recalling the Bianchi identity $\epsilon^{\alpha\beta\mu\nu}\partial_{\beta}F_{\mu\nu}=0$ . In this way, modifications to Maxwell equations only arise when $\partial_{\mu}\vartheta\neq 0$, as it happens at the interface $\Sigma$ of two materials having different constant values of $\vartheta$, for example. In this case, the action (2) can be integrated yielding a $2+1$ action at the boundaries corresponding to the Chern-Simons term. This means that at the boundaries of a 3D TI we have a quantum anomalous Hall (QAH) effect associated with each value of $\vartheta$, with Hall conductivity given by $\sigma_{{\rm H}}=\vartheta e^{2}/2\pi h$. In this way, the contribution to the total Hall conductivity from the interface $\Sigma$ between a TI and a regular insulator is

[TABLE]

since two values of $\vartheta$ differing by an integer multiple $m$ of $2\pi$ are equivalent.

The half integer contribution to $\sigma^{\Sigma}_{H}$ is a bulk property, which allows us to distinguish this case from that of a 2D surface gapped crystal having $\sigma_{H}=Ne^{2}/h$, with $N$ integer, thus showing that both conditions are not topologically equivalent. When dealing with a TR-invariant TI in a region with no boundaries, the number $m$ remains undetermined. The integer part of $\sigma^{\Sigma}_{H}$ becomes resolved only in the presence of a boundary between two TIs with different values of $\vartheta$, when TR symmetry is broken by gapping the interface. In this way we provide an adiabatic transition between those two topologically inequivalent insulators and the value of $m$ depends on the specific properties of such breaking. Such TR symmetry breaking is usually realized by an external magnetic field across the interface or by a magnetic doping of the surface.

The existence of TIs was predicted in Refs. Kane-Mele 1 ; Kane-Mele 2 ; Bernevig-Hughes-Zhang and their observation in two-dimensional HgTe/CdTe quantum wells was reported in Ref. Koenig et al . Then the authors in Refs. Fu-Kane-Mele ; Moore-Balents ; Roy discovered that the topological characterization of the quantum spin Hall insulator state has a natural generalization in three dimensions. Shortly afterward, this behavior was predicted in several real materials Fu-Kane which included Bi1-xSbx as well as strained HgTe and $\alpha$-Sn. Subsequently, the experimental discovery of the first three-dimensional TIs in Bi1-xSbx was reported in Ref. Hsieh et al . Later, a second generation of TIs, such as Bi2Se3, Bi2Te3, and Sb2Te3, was identified theoretically in Ref. Xia et al and experimentally discovered in Refs. Xia et al ; Zhang et al . This led to the detection of a huge variety of TIs ANDO and to their subsequent classification in a periodic table where different classes of these materials can be identified, distinguishing them by the presence or absence of time-reversal, particle-hole, and chiral symmetry TI11 . A new class of TIs, called axion insulators (AXIs), has been recently proposed as a new arena to probe topological phases. They have the same bulk MEP $\vartheta=\pi$ as 3D TIs, with gapped both bulk and surfaces, where the topological index is protected by inversion symmetry, instead of time-reversal symmetry. They are expected to show up in magnetically doped TI heterostructures with magnetization pointing inward and outward from the top and bottom interfaces of the TI DiXiao ; Mogi1 . Also, the possibility of having a different class of intrinsic axion insulators, which do not require magnetic doping, has been investigated in a pyrochlore lattice in Ref. Varnava .

When time reversal is broken at the interface between a 3D strong TI and a regular one, either by the application of a magnetic coating and/or by doping the TI with transition metal elements, the opening of the gap in the surface states induces several exotic phenomena that can be tested experimentally. Among them we find the (QAH) effect, the quantized magneto-optical effect, topological magnetoelectric (TM) effect, and the image magnetic monopole effect, all of which are a direct consequence of the change in the parameter $\vartheta$ between the two phases. These effects can be also produced in AXIs.

The QAH effect has already been experimentally observed in thin films of the magnetic chromium-doped topological insulator $({\rm Bi,Sb})_{2}{\rm Te}_{3}$ Chang 1 . In Ref. Kou the observation of the QAH effect in extremely thin films of the magnetic topological insulator (Cr0.12Bi0.26Sb0.62Sb0.62)2Te3 is reported. The characteristic behavior of the QAH effect has also been shown experimentally in thin films of the topological insulator Crx(Bi1-ySby)2-xTe3, which were grown on semi-insulating InP(111) substrates using molecular beam epitaxy methods Checkelsky . Also, Ref. Chang 2 demonstrates a high-precision confirmation of the QAH state in V-doped (Bi,Sb)2Te3 films, which is a strong ferromagnetic TI. The QAH effect has also been observed in AXIs Mogi . Employing terahertz time-domain spectroscopy, quantized magneto-optical effects have been observed by measuring the topological Faraday and Kerr rotations in QAH states on Crx(Bi0.26Sb0.74)2-xTe3 magnetic TI films. In this work the authors also report the observation of the QAH together with an experimental indication of the TM effect Okada . Quantized Faraday and Kerr rotations in magnetic fields higher than 5 T were observed in the 3D topological insulator Bi2Se3, providing evidence of the TM effect by an indirect measurement of the value $\vartheta=\pi$ Wu . Reference Dziom reported a quantized Faraday rotation in high external magnetic fields when linearly polarized terahertz radiation passes through the two surfaces of a strained HgTe 3D TI. This constitutes a direct consequence of the TM effect, thus confirming axion electrodynamics as the effective theory describing the response of 3D strong TIs.

Our following discussion of the reversed VC radiation parallels the analysis of VC radiation in standard references Schwinger ; Jackson ; Panofsky . We deal with radiation by further extending to the time-dependent case the Green’s function method already developed in Refs. UrrutiaMartinCambiaso3 ; UrrutiaMartinCambiaso to study the static electromagnetic response of TIs.

The paper is organized as follows. In Sec. II, we review the effective theory describing the electromagnetic response of magnetoelectrics, to be called $\vartheta$-ED. Section III contains the calculation of the time-dependent Green’s function (GF) for the electromagnetic potential $A_{\mu}$ arising from arbitrary sources in the presence of two semi-infinite magnetoelectric media with different values of the parameter $\vartheta$, each having the same permittivity $\epsilon$, and separated by a planar interface that encodes the manifestation of the ME effect. This section also reports the results of the far-field approximation of the GF, which is required to deal with radiation. We give an analytic expression of the GF in this limit, whose detailed calculations are presented in the Appendix A. In this way, we determine how the phase of the GF in standard ED Schwinger changes due to the presence of the magnetoelectric media, thus explaining the origin of reversed VC radiation. In Sec. IV, we consider a charged particle moving with constant velocity $v>c/\sqrt{\epsilon}$ in the direction perpendicular to the interface between the two magnetoelectrics. Our formalism yields analytic results for the physical quantities involved. The far-field expressions for ${\mathbf{E}}$ and ${\mathbf{B}}$ are calculated, together with the angular distribution of the total radiated energy per unit frequency $d^{2}E/d\omega d\Omega$. The resulting angular distribution indicates the presence of reversed VC radiation, which is the most important conclusion of our work. Its presence is further illustrated in polar plots showing the angular distribution of the reversed VC radiation together with those corresponding to the standard forward VC radiation. In Sec. V, we calculate the total radiated energy per unit frequency $dE/d\omega$, both in the forward and the backward directions, together with the number of photons per unit length radiated in the backward direction. The power per unit frequency radiated in the backward direction is also obtained. Section VI is devoted to some numerical estimations of our results considering the topological insulator ${\rm TlBiSe}_{2}$ as one of the semi-infinite media and a regular insulator with the same permittivity as the other. We have chosen this material because its electric and topological properties are known from previous references. We emphasize that we are not doing any ab initio or experimental calculation, either on this material or in any other sector of the manuscript. Let us insist again that our approach is based on the effective electromagnetic response of magnetoelectrics encoded in $\vartheta$-ED. From this point of view, we assume that any material we refer to has been fully characterized, meaning that their macroscopic parameters have been already determined. Section VII comprises a concluding summary of our results. In the following, we set $\hbar=c=1$, we denote by $\eta^{\mu}{}_{\nu}$ the Minkowski metric with signature $(+,-,-,-)$, and we take the convention $\varepsilon^{0123}=1$ for the Levi-Civita symbol.

II $\vartheta$ Electrodynamics

Let us consider two semi-infinite magnetoelectric media separated by a planar interface $\Sigma$ located at $z=0$, filling the regions $\mathcal{U}_{1}$ and $\;\mathcal{U}_{2}$ of the space, as shown in Fig. (1). Motivated by the results in Ref. UrrutiaMartin , whereby the effects of the $\vartheta$-term are substantially enhanced with respect to the optical contributions when both $\vartheta$ media have the same permittivity, we take $\epsilon_{1}=\epsilon_{2}=\epsilon$. In our case, we want to suppress the transition radiation, which occurs whenever a charge propagates across two media with different permittivities TRANSITION . This radiation would be present for all particle velocities and would interfere with the VC radiation we are interested in, making unnecessarily complicated its theoretical discussion, together with hindering its possible experimental verification. Additionally, we assume that the parameter $\vartheta$ is piecewise constant taking the values $\vartheta=\vartheta_{1}$ in the region $\mathcal{U}_{1}$ and $\vartheta=\vartheta_{2}$ in the region $\mathcal{U}_{2}$. This is expressed as

[TABLE]

where $H(z)$ is the Heaviside function with $H(z)=1,\;\;z\geq 0\;$, and $H(z)=0,\;\;z<0.$

In this restricted case we can take the action for the effective theory describing the electromagnetic response of this media as

[TABLE]

where $\varrho$ ($\mathbf{J}$) are charge (current) densities and $e^{2}=\alpha$ is the fine-structure constant.

Let us emphasize that we are describing each medium by a constant MEP $\vartheta$ in the bulk, which has the value $\vartheta=\pi$ ($\vartheta=0$) for a topological insulator (regular insulator). While there are boundary effects arising from the gapping of the interface, the MEP $\vartheta=\pi$ of the bulk of the 3D topological insulator has a gauge-invariant and topological origin.

The electromagnetic fields $\mathbf{E}$ and $\mathbf{B}$ are related with the electromagnetic potentials $\Phi$ and $\mathbf{A}$ in the standard form

[TABLE]

as a consequence of the homogeneous Maxwell equations

[TABLE]

It can be shown that the term ${\mathbf{E}}\cdot{\mathbf{B}}$ is a total derivative, which tells us that there are no modifications to the dynamics in the bulk. In this way, all the effects induced by ${\cal L}_{\vartheta}$ arise on the interface and will manifest themselves as a consequence of the boundary conditions there. Performing the variation of the action (5) gives the following set of modified Maxwell equations:

[TABLE]

where $\mathbf{\hat{u}}$ is the outward unit vector normal to the region ${\cal U}_{1}$ and

[TABLE]

In the case of a TI located in region $2$ of Fig. 1 (${\vartheta}_{2}=\pi$), in front of a regular insulator ($\vartheta=0$) in region $1$, we have

[TABLE]

with $m$ being an integer depending on the details of the TR symmetry breaking at the interface. For definiteness, we deal with media having $\mu=1$. The novelty in the above equations is that they introduce additional field-dependent effective charge and current densities

[TABLE]

with support only on the interface $\Sigma$ between the two media. Consequently, the standard Maxwell equations remain valid in the bulk. The conservation equation

[TABLE]

can be readily verified by using Faraday’s law together with the relation

[TABLE]

which arises since the only nonzero component of $u_{k}$ is $u_{3}=1$.

The terms proportional to $\tilde{\theta}$ in Eqs. (8) and (9) describe the ME effect which is the distinctive feature of $\vartheta$-ED. Let us remark that Eqs. (8) and (9) can also be obtained starting from the standard Maxwell equation in a material medium Schwinger ; Jackson ,

[TABLE]

with the constitutive relations

[TABLE]

Assuming that the time derivatives of the fields are finite in the vicinity of the interface, the modified Maxwell equations (8) and (9) yield the following boundary conditions (BCs)

[TABLE]

for vanishing external sources at $z=0$. These BCs are derived either by integrating the field equations over a pill-shaped region across the interface or by using the Stokes theorem over a closed rectangular circuit perpendicular to the interface Schwinger ; Jackson . The notation is $\left[\mathbf{V}\right]_{z=0^{-}}^{z=0^{+}}=\mathbf{V}(z=0^{+})-\mathbf{V}(z=0^{-})$, \mathbf{V}\big{|}_{z=0}=\mathbf{V}(z=0), where $z=0^{\pm}$ indicates the limits $z=0\pm\eta$, with$\;\eta\rightarrow 0$,  respectively. Here the subindex $\parallel$ denotes the projection of a vector into the $x-y$ plane. The continuity conditions (19) imply that the right-hand sides of the discontinuity conditions (II) are well defined and they represent self-induced surface charge and surface current densities, respectively. The BCs (II) clearly illuminate again the ME effect, which is localized just at the interface $\Sigma$ between the two media.

III Green’s function method

In this section, we extend to the case of time-dependent $\vartheta$-ED with planar symmetry the Green’s function method discussed in Ref. UrrutiaMartinCambiaso for the static limit. A first step in this direction was presented in Ref.UrrutiaMartinCambiaso3 . The knowledge of the GF allows us to compute the electromagnetic fields for an arbitrary distribution of sources, as well as to solve problems with given Dirichlet, Neumann, or Robin boundary conditions on arbitrary surfaces. In what follows, we restrict ourselves to contributions of free external sources $J^{\mu}=(\varrho,\mathbf{J})$ located outside the interface $\Sigma$ and to systems without boundary conditions imposed on additional surfaces, except for the boundary conditions at infinity.

A compact formulation of the problem can be given in terms of the potential $A^{\mu}=(\Phi,\mathbf{A})$, which satisfies the equations

[TABLE]

in the modified Lorenz gauge

[TABLE]

with the operator

[TABLE]

The current is $J^{\mu}=(\varrho,{\mathbf{J}})$. The BCs (II) and (19) reduce to

[TABLE]

in terms of the vector potential. Next we introduce the Green’s function $G_{\;\;\sigma}^{\nu}(x,x^{\prime})$ satisfying

[TABLE]

together with the BCs arising from Eq. (23), in such a way that the four potential is

[TABLE]

determined up to homogeneous solutions of Eq. (20). We solve Eq. (24) along the same lines introduced in Refs. UrrutiaMartinCambiaso ; UrrutiaMartinCambiaso3 . First we take advantage of the translational invariance in time and in the transverse directions $x$ and $y$ by introducing the reduced Green’s function $g_{\;\;\nu}^{\mu}(z,z^{\prime};\mathbf{k_{\perp}},\omega)$, such that Schwinger

[TABLE]

where $\mathbf{R}_{\perp}=(\mathbf{x}-\mathbf{x^{\prime}})_{\perp}=(x-x^{\prime},y-y^{\prime})$ and $\mathbf{k_{\perp}}=(k_{x},k_{y})$ is the momentum perpendicular to the vector ${\mathbf{\hat{u}}}$ in Fig. 1. The resulting equation for $g_{\;\;\nu}^{\mu}$ is

[TABLE]

where  $k^{\alpha}=(\omega,\mathbf{k}_{\perp},0)$, and ${\mathcal{O}}^{\mu}{}_{\nu}$ is given by Eq. (22) after the replacements $\nabla^{2}\rightarrow-{\mathbf{k}}^{2}_{\perp}+\partial^{2}/\partial z^{2}$ and $\partial_{t}\rightarrow-i\omega$ are made.

To solve Eq. (27) we employ the — by now — standard method Grosche to deal with $\delta$-like interactions, whereby a convenient free GF can be used to integrate the GF equation (27). By free, here we mean a GF satisfying Eq. (24) with $\tilde{\theta}=0$ and which can be written in the same form as in Eq. (26). To this end, we consider the reduced free Green’s function ${\mathcal{H}}^{\rho}{}_{\nu}$ associated with the operator ${\mathcal{O}}^{\mu}{}_{\nu}$ previously defined, which solves

[TABLE]

and satisfies standard boundary conditions at infinity. From now on, calligraphic letters will denote free GFs, i.e., GFs that are independent of $\vartheta$. Separating the components, we have

[TABLE]

which implies

[TABLE]

and leaves us with only

[TABLE]

This system can be solved in terms of the function ${\mathcal{F}_{0}}(z,z^{\prime};\mathbf{k_{\perp}},\omega)$ satisfying

[TABLE]

plus the standard BCs at infinity. The result is Schwinger

[TABLE]

with

[TABLE]

where we define $\tilde{k}_{0}=\omega\sqrt{\epsilon}$ with $\omega>0$ and $\tilde{k}^{\alpha}=(\sqrt{\epsilon}\omega,\mathbf{k}_{\perp},0),\;\tilde{k}^{\alpha}\tilde{k}_{\alpha}=\epsilon\omega^{2}-\mathbf{k}_{\perp}^{2}=\tilde{k}^{2}$. In this way, we obtain

[TABLE]

Then, Eq. (27) can be integrated using ${\mathcal{H}}_{\;\;\sigma}^{\nu}(z,z^{\prime})$ as follows

[TABLE]

yielding

[TABLE]

An example of the flexibility of this method is given in Ref. UrrutiaMartinCambiaso3 , where the free GF was chosen to satisfy the required BCs in order to limit the space in direction $z$ to the region between two metallic plates parallel to the interface. Equation (39) reduces the calculation of $g_{\;\;\sigma}^{\nu}$ to solving a set of linearly coupled algebraic equations. The details of the procedure for solving Eq.(39) are presented in Ref. UrrutiaMartinCambiaso . The solution is

[TABLE]

where

[TABLE]

and

[TABLE]

with $u^{\mu}=(0,0,0,1)$. It is interesting to observe that the term $P(z,z^{\prime};\mathbf{k_{\perp}},\omega)$ is proportional to $e^{ik_{z}(|z|+|z^{\prime}|)}$, whose phase is finally responsible for the reversed VC radiation, as we will show in Sec. IV.

Since the only difference between the GFs $G^{\mu}{}_{\nu}$ and ${\bar{G}}^{\mu}{}_{\nu}$ is that $G^{0}{}_{0}={\bar{G}}^{0}{}_{0}/\epsilon$, with all other terms being equal, it is convenient to present the barred GFs that are obtained once $\bar{g}_{\;\;\nu}^{\mu}(z,z^{\prime};\mathbf{k_{\perp}},\omega)$ is substituted in Eq. (26). The full GF matrix $\bar{G}_{\;\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)$ can finally be written as the sum of three terms, $\bar{G}_{\;\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)=\bar{G}_{ED\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)+\bar{G}_{\tilde{\theta}\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)+\bar{G}_{\tilde{\theta}^{2}\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)$, where the three pieces are

[TABLE]

Let us emphasize that the potentials and the resulting electromagnetic fields must be calculated using the GF ${G}_{\;\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)$. In the static limit ($\omega=0$), the result (43) reduces to that reported in Refs. UrrutiaMartinCambiaso .

The next step is to evaluate the GFs (43) in the far-field approximation corresponding to the coordinate conditions

[TABLE]

where $\|\mathbf{x}\|=r\rightarrow\infty$, $\rho\rightarrow\infty$, and $z\rightarrow\infty$. In this approximation, the integrals in Eqs. (43) include rapidly oscillating exponential functions whose leading contributions are calculated in the stationary phase approximation and subsequently verified by the steepest descent method Chew ; Chew1 ; Mandel . Also we make repeated use of the generic Sommerfeld identity Sommerfeld

[TABLE]

where ${\cal R}=\sqrt{R_{\perp}^{2}+Z^{2}}$. We consider mainly two cases: (i) $Z=|z-z^{\prime}|$ where ${\cal R}$ is denoted by $R$ and (ii) $Z=|z|+|z^{\prime}|$ where ${\cal R}$ is denoted by ${\tilde{R}}$. The detailed calculation is presented in the Appendix and the results in the far-field approximation are

[TABLE]

where we define

[TABLE]

Here ${s}_{\alpha}=(1/n,\hat{\mathbf{n}})$, $\hat{\mathbf{n}}$ is a unit vector in the direction of $\mathbf{x}$, and $\mathbf{n}_{\perp}=(x/r,\,y/r\,,0)=\sin\theta(\cos\phi,\,\sin\phi,\,0)$ is a vector in the direction of $\mathbf{x}_{\perp}$ with $\|\mathbf{n}_{\perp}\|=\sin\theta$. In Eq. (49) we have introduced the shorthand notation $x=r\sin\theta\cos\phi$ and $y=r\sin\theta\sin\phi$ in spherical coordinates.

As explained in the Appendix, the results for the far-field approximation of the GFs in Eqs. (46)–(48) are valid whenever $0<\theta<{\theta}_{\rm C}=\cos^{-1}(1/vn)<\cos^{-1}(1/n)\ll\pi/2$, in such a way that there are no divergences arising from the factors proportional to $1/\cos\theta$ in the GFs, since $\theta$ is never equal to $\pi/2$.

The GFs (46)–(48) together with Eqs. (6) and (25) yield electromagnetic fields whose Cartesian components behave like $(e^{i{\tilde{k}}_{0}r}/r)F(\theta)$ in the far-field approximation. Recalling that ${\tilde{k}}_{0}=\omega\sqrt{\epsilon}$ and using Maxwell equations to leading order in $1/r$ we verify the expressions

[TABLE]

which are the distinctive feature of the radiation fields. The three vectors ${\mathbf{E}}$, ${\mathbf{B}}$, and ${\mathbf{\hat{n}}}$ define a right-handed triad resulting in the Poynting vector for a material media with $\mu=1$,

[TABLE]

where ${\mathbf{\hat{n}}}$ is in the direction of the phase velocity of the outgoing wave, as appropriate for right-handed materials.

Finally, it is pertinent to emphasize an important difference in the phase of the exponentials related to the source variables ${\mathbf{x}}^{\prime}$ in the GFs (46)–(48). In the first case, we encounter the exponential $e^{i{\tilde{k}}_{0}R}$, which in the coordinate approximation of the far-field zone produces the phase $i{\tilde{k}}_{0}(r-{{\mathbf{n}}}_{\perp}\cdot{\mathbf{x}}^{\prime}_{\perp}-z^{\prime}\cos\theta)$ characteristic of radiation in standard electrodynamics Schwinger . On the other hand, the contributions to the GF proportional to $\tilde{\theta}$ and $\tilde{\theta}^{2}$ involve the exponential $e^{i{\tilde{k}}_{0}{\tilde{R}}}$ with ${\tilde{R}}=\sqrt{({\mathbf{x}}-{\mathbf{x}}^{\prime})_{\perp}^{2}+(|z|+|z^{\prime}|)^{2}}$, which can be presented as follows:

[TABLE]

in the far-field approximation. From the second line in Eq. (52), we remark that whenever the sign of $zz^{\prime}$ is positive we will have an additional relative phase contributing to the GFs, which will show up in observable quantities as the radiated power, for example. The term $(|z|+|z^{\prime}|)^{2}$ can ultimately be traced back to the form of the reduced GF (41) together with expressions (34) and (42). As we will show in the next section, reversed VC radiation arises precisely due to the contribution $|z^{\prime}\cos\theta|$ in the phase of the GF deriving from the term $(|z|+|z^{\prime}|)^{2}$.

IV The reversed Vavilov-Čerenkov radiation

Let us now consider a particle with charge $q$ moving with constant velocity $v{\mathbf{\hat{u}}}$, perpendicular to the interface $\Sigma$ defined by the $x-y$ plane ($z=0$), as shown in Fig.1. The charge and current densities are

[TABLE]

where we henceforth assume $v>1/\sqrt{\epsilon}>0$. Instead of an infinite path for the charge, we will take its movement in the interval $z\in(-\zeta,\zeta)$, with $\zeta\gg v/\omega$. In the far-field approximation, the resulting components of the electric field $\mathbf{E}(\mathbf{x};\omega)=-i{\tilde{k}}_{0}\mathbf{\hat{n}}A^{0}(\mathbf{x};\omega)+i\omega\mathbf{A}(\mathbf{x};\omega)$, calculated in terms of the potential $A^{\mu}(\mathbf{x};\omega)$ via Eq. (25), are

[TABLE]

We can verify that $\mathbf{{\hat{n}}\cdot E}=0$ as required. The magnetic field can be obtained from Eq. (50). The integrals $\mathcal{I}_{1}(\omega,\theta)$ and $\mathcal{I}_{2}(\omega,\theta)$ appear when taking the convolution of the charge and current densities with the GF components in Eqs. (46)–(48). They are defined by

[TABLE]

where

[TABLE]

In the following, we will learn that the production of VC radiation depends on the zeros of $\Xi_{\pm}$ and ${\tilde{\Xi}}_{\pm}$, which determine the angle $\theta_{C}$ of the VC cone. Then it is reasonable to expect that the $|\cos\theta|$ dependence of ${\tilde{\Xi}}_{\pm}$ yields new possibilities. In calculating the right-hand side of Eqs. (57) and (58) in the limit $\zeta\gg v/\omega$, which effectively means $\zeta\rightarrow\infty$, we encounter expressions like $\sin(\zeta aN)/(aN)$, which behave as $\pi\delta(aN)$ Lighthill . We take advantage of this $\delta$-like behavior by setting equal to zero all the rapidly oscillating contributions arising from those functions with an argument that can never be zero, like $\tilde{\Xi}_{+}$, for example. This is relevant to obtain the final expression for $\mathcal{I}_{2}(\omega,\theta)$, which then simplifies to

[TABLE]

Recalling that $\mathbf{B}=\hat{\mathbf{n}}\times\left(n\mathbf{E}\right)$ for radiation fields, with $n=\sqrt{\epsilon}$, the angular distribution of the total radiated energy per unit frequency in the interval $-\zeta<z<+\zeta$ with $\zeta\rightarrow\infty$ is Schwinger ; Jackson

[TABLE]

which can be written as the sum of the following three terms:

[TABLE]

To obtain the above expressions we have used again the $\delta$-like behavior of the functions ${\cal I}_{1}$ and ${\cal I}_{2}$ by replacing $\theta$ by $\theta_{C}=\cos^{-1}(1/vn)$ in all the functions of $\theta$ that multiply them. We have verified the cancellation of the terms proportional to ${\tilde{\theta}}^{3}$. Also, the addition of the contributions in the term containing the factor $1/(4n^{2}+{\tilde{\theta}}^{2})^{2}$ is such that the final result ends up being proportional to ${\tilde{\theta}}^{2}/(4n^{2}+{\tilde{\theta}^{2})}$.

The quantity $d^{2}E/d\omega d\Omega$ is also referred to as the spectral distribution of the radiation Schwinger , a shorter synonymous that we will use in the following. Let us observe that the above distributions have azimuthal symmetry and are even functions of both the angle $\theta$ [recall Eqs. (57) and (60)] and the MEP $\tilde{\theta}$. In other words, the leading corrections arising from the ME effect are of order ${\tilde{\theta}}^{2}$. Setting $\tilde{\theta}=0$ we recover the well-known expression for the spectral distribution of the radiation in the standard VC case Panofsky

[TABLE]

Substituting $\mathcal{I}_{1}(\omega,\theta)$ and $\mathcal{I}_{2}(\omega,\theta)$ from Eqs. (57) and (60) in Eqs. (63) and (LABEL:d2E_2), we have

[TABLE]

Equations (65)–(67) summarize the spectral distribution of the radiation when the charges passes through a $\vartheta$ medium. The $\delta$-like behavior of the functions appearing there yields the VC radiation condition $\sin^{2}\theta_{C}=1-1/(v^{2}n^{2})$ and clearly shows that the distributions (65) and (66) contribute only to the forward Vavilov-Čerenkov radiation with $\cos\theta_{C}=1/(nv)>0$. This is precisely due to the dependence upon $|\cos\theta|$ in the angular distribution originating from ${\tilde{\Xi}}_{-}$ in Eq. (59).

Even though the spectral distribution in Eq. (66) depends on ${\tilde{\Xi}}_{-}$, its contribution shows up in a product that behaves like $\delta(1-vn\cos\theta)\times\delta(1-vn|\cos\theta|)$ in the limit of large $\zeta$. This yields a nonzero spectral distribution $d^{2}E_{12}/d\omega d\Omega$ only when $\cos\theta>0$, i.e., in the forward direction. On the contrary, the distribution $d^{2}E_{2}/d\omega d\Omega$ contributes both to the forward VC cone ($\cos\theta>0$), as well as to the backward VC cone ($\cos\theta<0$) because ${\tilde{\Xi}}_{-}$, depending upon $|\cos\theta|$, admits also a zero in the range $\pi/2<\theta<\pi$. In other words, radiation is also detected in the backward direction, i.e., $\cos\theta<0$, according to the angular distribution (67). This defines the reversed VC cone and constitutes the most important result of the manuscript.

The forward VC radiation receives corrections of order ${\tilde{\theta}}^{2}$ with respect to the standard case. The reversed VC radiation is purely of order ${\tilde{\theta}}^{2}$ and, in general, will be strongly suppressed with respect to the forward output, however, it is different from zero. If the angle of the forward cone is $\theta_{0}$, the one corresponding to the reversed cone will be $\pi-\theta_{0}$. Since this is standard electrodynamics plus additional field-dependent sources at the interface, the group and phase velocities of the electromagnetic wave are parallel, as already pointed out after Eq.(51). In fact this model provides the usual interpretation of the radiation of a system observed at an arbitrary solid angle $d\Omega$ Schwinger ; Jackson .

Figures 2 and 3 show the spectral distribution in the case of pure forward VC radiation, Eq. (65), in a dielectric with $n=2$ and $\tilde{\theta}=0$, for the average frequency of $\omega=2.5$ eV ($500$ nm) in the VC radiation spectrum. In each figure, the particle moves from left to right along the line ($\pi-0$). In Fig. 2, the dashed purple line corresponds to $v=0.9$ and $\zeta=1.0$ eV*-1* and the solid cyan line to $v=0.5009$ and $\zeta=2.8$ eV*-1*. Such low values of $\zeta$ are chosen to illustrate the emergence of the forward radiation lobes. In Fig. 3, the dashed purple line corresponds to $v=0.9$ and $\zeta=340$ eV*-1* and the solid cyan line to $v=0.5009$ and $\zeta=4800$ eV*-1*. Such an increase in each value of $\zeta$ is enough to show the appearance of the forward VC cone in both cases.

Figures 4 and 5 include the spectral distribution of the radiation arising from the contribution of Eq. (67) in a medium with ${\tilde{\theta}}\neq 0$, which clearly show the presence of reversed VC radiation . This term also provides corrections to the forward VC radiation. Both additions are highly suppressed with respect to the forward VC radiation, so that the scales in both figures are separately chosen such as to make these small, but nonzero, contributions clearly visible. Here we take again the frequency of $\omega=2.5$ eV and the charge propagates now in a $\vartheta$ medium with $n=2$ and $\tilde{\theta}=11\alpha$, from left to right along the line ($\pi-0$). In Fig. 4, the dashed purple line corresponds again to $v=0.9$ and $\zeta=1.0$ eV*-1* and the solid cyan line to $v=0.5009$ and $\zeta=2.8$ eV*-1*. The contribution of $d^{2}E_{2}/d\omega d\Omega$ to the forward and reversed VC cones are displayed in Fig. 5, where the parameters $\zeta$ have been modified with respect to those in Fig. 4 by increasing them to $343$ eV*-1* (dashed purple line) and $4830$ eV*-1* (solid cyan line), respectively. All the spectral distributions we have plotted are calculated from the respective expressions $d^{2}E/d\omega d\Omega$ and the results are expressed in units of the common factor $q^{2}/\pi^{2}=7.4\times 10^{-4}$ for $q=\sqrt{\alpha}$. Let us emphasize that the negative contribution in Eq. (66) only diminishes the radiation in the forward direction but does not affect the reversed VC radiation.

The full spectral distribution of the total radiated energy is given by the sum of Eqs. (65)–(67) and it is plotted in Fig. 6. The scale of the polar plot is in arbitrary dimensions and runs from [math] to $10^{6}$. On the left side of the figure we plot an elargement in the backward direction showing the onset of the reversed VC radiation. Here the radial scale goes from 0 to $10^{2}$ showing the high suppression of the radiation in the backward direction with respect to that in the forward direction.

V Total radiated energy

In this section, we calculate the total energy per unit frequency radiated by the charge on its path from $-\zeta$ to $+\zeta$. Let us first review the calculation of ${dE_{1}}/{d\omega}$, following the procedure in Ref. Panofsky . Integrating expression (65) with respect to the solid angle we obtain

[TABLE]

The $\delta$-like behavior of the integrand in the limit $\zeta\gg\omega/v$ shows that the radiation is sharply localized in a main lobe around the angle $\theta_{C}$ given by $\cos\theta_{C}=1/(vn)$, yielding a contribution only in the forward direction. Therefore, making the change of variable $u=\frac{\omega\zeta}{v}\left(1-vn\cos\theta\right)$, we can safely extend the integration limits of $u$ to $\pm\infty$ as long as we include the maximum located at $u=0$. The result is

[TABLE]

where we have introduced the total length $L=2\zeta$ traveled by the particle, thus recovering the standard result Panofsky . In other words, we are estimating the contribution from each sharply localized lobe as

[TABLE]

The next term comes from the angular integration of Eq. (66),

[TABLE]

The presence of $|\cos\theta|$ from the function $(\zeta\tilde{\Xi}_{-})$ in the integrand requires separating the integration with respect to $\theta$ into the regions $0<\theta<\pi/2\,\,$ and $\,\,\pi/2<\theta<\pi$. After splitting the integral as required and setting equal to zero the rapidly oscillating contribution proportional to $\sin\left(\zeta\Xi_{-}\right)$ in the region $\pi/2<\theta<\pi$, we are left with

[TABLE]

Again, this contribution to the radiation is only in the forward direction and it is concentrated in a main lobe around $\theta_{C}$. Making use of Eq. (70) yields

[TABLE]

The calculation of ${dE_{2}}/{d\omega}$ starts from Eq. (67), whose separate contributions are

[TABLE]

The presence of $|\cos\theta|$ in $\tilde{\Xi}_{-}$ indicates that now the radiation is concentrated in two main lobes, one around $\theta_{C}$, contributing to the forward radiation, and the other around $\pi-\theta_{C}$, contributing to the reversed VC radiation. Taking into account both lobes when integrating Eq.(74) over the solid angle and using Eq. (70), we find

[TABLE]

The term ${d^{2}E_{2b}}/{d\omega d\Omega}$ in Eq. (75) also contributes both to the forward and backward radiation and differs from the case of (74) by the replacements $\tilde{\Xi}_{-}\rightarrow\tilde{\Xi}_{-}/2$ and $\sin^{2}u\rightarrow\sin^{4}u$. Nevertheless, we obtain the same result as in Eq. (76),

[TABLE]

Therefore, the total radiated energy per unit frequency is

[TABLE]

which corresponds to the same expression as in the absence of the $\vartheta$ medium. At this stage, we do not know whether Eq. (78) is just a coincidence or there is a fundamental reason for this result. For our purposes, the main point is that the total energy distribution (78) is split into a nonzero reversed VC radiation (RVCR)

[TABLE]

together with a modified forward VC radiation (FVCR)

[TABLE]

We can restate expression (79) in terms of the number of photons radiated per unit length and per unit frequency as

[TABLE]

Another relevant parameter is the power radiated per unit frequency in the backward direction,

[TABLE]

whose estimation in Table 1 will be used for further comparison with alternative predictions.

VI Order of magnitude evaluations

In order to present some numerical estimations we use the setup in Fig. 1 where medium 1 is a regular insulator with $\epsilon=4,\mu=1,\vartheta=0$, and medium 2 is the topological insulator ${\rm TlBiSe}_{2}$ with $\epsilon=4,\,n=2,\,\mu=1$, ${\vartheta}=\pi$, and $0<m<5$ TlBiSe2 . This choice of $m$ yields a value of $\tilde{\theta}$ in the range $[\alpha,11\alpha]$, which very much suppresses the backward radiation as compared with the forward output. Our numerical estimations do not consider the frequency dependence of the refraction index. In this way, we can calculate the ratio

[TABLE]

where we also take $4n^{2}\gg{\tilde{\theta}}^{2}$. For the topological insulator ${\rm TlBiSe}_{2}$ the above ratio ranges between $1.7\times 10^{-6}$ and $2.0\times 10^{-4}$ for the choices ${\tilde{\theta}}=\alpha$ and ${\tilde{\theta}}=11\alpha$, respectively. In Table 1, we show the power radiated per unit frequency in the backward direction, given by Eq. (82), for ${\rm TlBiSe}_{2}$ with $\tilde{\theta}=11\alpha$. The results are in units of $\left[{\mu W}/{eV}\right]$, for $\omega\in[2,\,8]$ eV. This range includes the major sector of the frequency spectrum in the forward VC radiation. The conversion to MKS units is 1$\,\mu$W/eV= 6.24$\times$ 1012 s*-1*, which is equal to $4.11\times 10^{-4}$ eV in units where $c=\hbar=1$.

VII Summary

We have considered the radiation produced by an electric charge propagating with constant velocity $v$ in the direction ${\mathbf{\hat{u}}}$ (shown in Fig. 1) between two $\vartheta$ media with the same permittivity, whose electromagnetic response is driven by the modified Maxwell equations (7)–(10). When $v$ is higher than the speed of light in the media, we discover the emission of reversed VC radiation, codified in the angular distribution given by Eq. (67) and illustrated in Figs. 4–6. These right-handed $\vartheta$ media are realized in nature as magnetoelectric materials, among which we find topological insulators, having positive permittivity, permeability, and index of refraction.

The main characteristics of the reversed VC radiation we have discovered are the following: (i) The threshold condition $v>c/n$ for the velocity of the charge must be satisfied as in the standard case. (ii) The reversed VC radiation occurs for all frequencies in the VC spectrum and it is always accompanied by forward VC radiation. (iii) The energy loss per unit frequency of the reversed VC radiation is highly suppressed with respect to the forward output according to Eq. (83). A comparison with measurements of reversed VC radiation in metamaterials can be established by interpreting this suppression as due to the detection of radiation at an effective frequency $\omega_{\rm eff}=\omega\,{\tilde{\theta}}^{2}/8n^{2}$, according to Eq. (79). Standard Čerenkov counters work in the range of $140-800\,{\rm nm}$ corresponding to $8.9-1.6\,\,{\rm eV}$, respectively. Taking an average of $500\,\,{\rm nm}$ ($2.5\,\,{\rm eV}$), we would expect detectable reversed VC radiation at $\omega_{\rm eff}$ in the range from $4\times 10^{-3}\,{\rm meV}$ for ${\tilde{\theta}}=\alpha$ to $0.5\;{\rm meV}$ for ${\tilde{\theta}}=11\alpha$, respectively, using ${\rm TlBiSe}_{2}$ as a $\vartheta$ medium. However, recent measurements of reversed VC radiation in metamaterials show that these estimations are within the experimental capabilities. Reversed VC radiation has been measured at a frequency of $2.85\,\,{\rm GHz}$, equivalent to $1.2\times 10^{-2}\,\,{\rm meV}$, in an all-metal metamaterial consisting of a square waveguide loaded with complementary electric split-ring resonators Duan . Likewise, reversed VC radiation in the range $(3.4-3.9)\times 10^{-2}\,\,{\rm meV}$ has also been experimentally verified in a phased electromagnetic dipole array used to model a moving charged particle Sheng .

Our estimations for $d^{2}E_{\rm RVCR}/dtd\omega$ in Table 1 are in the range $10^{-3}-10^{-2}\,\,\mu{\rm W/eV}$, in the frequency interval of $2-8\,\,{\rm eV}$. They are smaller by a factor of $10^{-4}-10^{-3}$ than the maximum output of $\sim 10\,\mu{\rm W/eV}$ theoretically predicted to occur in the narrow interval of $5.7-6.5$ eV in a metal-insulator-metal waveguide Tao . In such a waveguide with a core thickness of $a=20\,\,{\rm nm}$, surface plasmon polaritons excited by an electron moving at $v=0.8$ produce reversed VC radiation.

A qualitative argument for the existence of the reversed Čerenkov radiation in a $\vartheta$ media can be given by extending the interpretation of the static fields produced by a point charge located in front of a topological insulator in terms of electric and magnetic images UrrutiaMartinCambiaso ; QI_SCIENCE . When the charge $q$ moves from region 1 to region 2, as we have assumed here, the effect of the interface $\Sigma$ can be replaced by introducing a moving image electric charge ${\tilde{q}}=-q{\tilde{\theta}}^{2}/(4n^{2}+{\tilde{\theta}}^{2})$ together with a moving image magnetic monopole ${\tilde{g}}=2q{\tilde{\theta}}/(4n^{2}+{\tilde{\theta}}^{2})$, both located in region 2 and moving towards region 1 . These images would contribute to the physical fields only in region 1 with their own forward VC radiation in the far-field zone, which turns out to be in the reversed direction with respect to the incident charge. As shown in Ref. Schwinger , a magnetic monopole $g$ propagating with a velocity $v>1/n$ also produces a forward VC cone, with radiated energy per unit length and per unit frequency given by

[TABLE]

in complete analogy with Eq. (69). From this point of view, the dominant contribution to the reversed VC radiation of the incident charge would arise from the image magnetic monopole since $g\sim{\tilde{\theta}}$ and the image electric charge would contribute with higher order terms of order ${\tilde{q}}^{2}\sim{\tilde{\theta}}^{4}$. A detailed calculation of the radiation fields in region 1 produced by this configuration is further required to test this interpretation, in particular, to verify that the final factor depending on $\tilde{\theta}$ will have the correct form given in Eq. (79). As it is well known, image charges are only useful mathematical tools in the description of electromagnetic phenomena and do not represent physically existing entities. In our case, the physical response of the medium is produced by the time-varying electric charge densities and Hall currents induced at the interface $\Sigma$ due to the magnetoelectric effect. The verification of the above interpretation of the reversed VC radiation is beyond the scope of the present paper and is postponed for future work.

Our results apply for any material that supports the magnetoelectric effect, in such a way that its electromagnetic response can be described by $\vartheta$-ED. In the particular case of a topological insulator, the conditions for the realization of such effect are (i) the topological insulator should be in the 3D regime, (ii) all the surfaces need to be gapped with the chemical potential lying within the gaps and (iii) the dynamics of the interior of the topological insulator should be invariant under time-reversal symmetry or inversion symmetry, in order to keep $\vartheta=\pi$ in the bulk DiXiao . In this work, we have highlighted the case of 3D strong topological insulators, which require time-reversal symmetry breaking in the surfaces to realize the topological magnetoelectric effect. This can be achieved by magnetically gapping all the interfaces, in such a way that the entire sample behaves as an insulator having a magnetoelectric coupling of exactly $\vartheta=\pi$ VDB5 . Since we have assumed that charged particles are moving with constant velocity across the material, it would be advisable to break the time-reversal symmetry by doping the surfaces with thin ferromagnetic films instead of switching on an external magnetic field. In any case, the velocity of such moving particles has the lower limit of $c/n$ but must be large enough so that there is no appreciable deflection of the charge when going through the magnetically doped material. The velocity of the external moving charges is independent of the Fermi velocity of the electrons on the 2D surfaces and the Fermi energy can comfortably lie in the middle of the gap between the Dirac cones. We expect our results to hold also in the case of magnetically doped AXIs. They are heterostructures in which magnetic ions are added to the vicinity of the top and bottom interfaces of a 3D strong topological insulator, like (Bi,Sb)2Te3, for example. In this way, their electromagnetic response is still coded in $\vartheta$-ED. Even though the upward-downward magnetic coating in opposite interfaces will produce a change in the sign of $\vartheta=\pi$ when going from one interface to the other, this will make the contribution of both interfaces to the spectral distribution of the reversed VC radiation to add up. This is because each contribution depends on ${\tilde{\theta}}^{2}$, respectively, according to Eqs. (11) and (67).

Finally, we comment on the condition $\epsilon_{1}=\epsilon_{2}$, imposed at the beginning of Sec. II in order to get rid of the otherwise unavoidable transition radiation and aimed to present a clean derivation of the reversed VC radiation. We do not mean that relaxing the equal permittivity condition will eliminate the presence of reversed VC radiation, but we only claim that the choice $\epsilon_{1}$ very different from $\epsilon_{2}$ will make transition radiation interfere with the VC radiation requiring a new theoretical discussion as well as hindering any experimental detection. In fact, the spectral distribution of the transition radiation might also display unexpected corrections arising from $\vartheta$-ED. In realistic terms, the equal permittivity condition means that $\epsilon_{1}$ and $\epsilon_{2}$ could be chosen such that their values are as close as possible. This condition can be considered as a useful restriction in the planning of a possible experimental setup and looks plausible since the refraction indices of normal insulators cover a wide range of values.

Acknowledgements.

O. J. F. has been supported by the doctoral fellowship CONACYT-271523. O. J. F., L. F. U. and O. R. T. acknowledge support from the CONACYT Project No. 237503. Support from the project No. IN103319 from Dirección General de Asuntos del Personal Académico (Universidad Nacional Autónoma de México) is also ackowledged. L. F. U. and O. J. F. thank Dr. Alberto Martín-Ruiz, Professor Hugo Morales-Técotl and Professor Rubén Barrera for useful discussions and suggestions. The authors also thank Dr. Christine Gruber for a careful reading of the manuscript together with many helpful comments.

Appendix A Green’s function in the radiation zone

In this Appendix, we obtain the far-field approximation for the GFs in Eq. (43) through the stationary phase method Chew . We review the contribution of standard ED and consider only ${\bar{G}}^{\mu}_{\tilde{\theta}\;\nu}(\mathbf{x},\mathbf{x}^{\prime};\omega)$ to illustrate the procedure of how to deal with the additional contributions arising from ${\cal L}_{\vartheta}$.

Let us start with the components ${\bar{G}}_{ED\;\nu}^{\mu}(\mathbf{x},\mathbf{x}^{\prime};\omega)$. The double integral in $\mathbf{k}_{\perp}$ is conveniently calculated by expressing the area element in polar coordinates $d^{2}\mathbf{k_{\perp}}=k_{\perp}dk_{\perp}d\varphi$ and choosing the $k_{\perp x}$ axis in the direction of the vector $\mathbf{R}_{\perp}=(\mathbf{x}-\mathbf{x^{\prime}})_{\perp}$. This defines the coordinate system ${\cal S}$ to be repeatedly used in the following and shown in the Fig. 7. Writing $\mathbf{k_{\perp}}\cdot\mathbf{R}_{\perp}=k_{\perp}R_{\perp}\cos\varphi$ and recalling that the angular integral provides a representation of the Bessel function $J_{0}(k_{\perp}R_{\perp})$ Abramowitz we obtain

[TABLE]

with $R=\|\mathbf{x}-\mathbf{x}^{\prime}\|=\sqrt{R_{\perp}^{2}+(z-z^{\prime})^{2}}$, where the right-hand side of the above equation is a direct consequence of the Sommerfeld identity Sommerfeld . Moreover, the coordinate conditions

[TABLE]

in the far-field approximation produce the well-known result for ${\bar{G}}^{\mu}_{ED\;\nu}$ in standard ED Schwinger ,

[TABLE]

with $\hat{\mathbf{n}}$ being a unit vector in the direction of $\mathbf{x}$ and $\|\mathbf{x}\|=r$ .

For future purposes, it will be convenient to go back to the generalized Sommerfeld identity in Eq.(45) and to rewrite it in terms of Hankel functions. We start from

[TABLE]

where $H_{0}^{(1)}$ and $H_{0}^{(2)}$ are the Hankel functions, together with the reflection formula $H_{0}^{(1)}(e^{i\pi}x)=-H_{0}^{(2)}(x)$, which allows us to extend the integration interval to $-\infty$. The result is

[TABLE]

with ${\cal R}=\sqrt{R_{\perp}^{2}+Z^{2}}$ and where $Z$ can be conveniently chosen. Here $C$ denotes the path of integration defined in Fig. 2.2.5 of Ref. Chew .

Now, we calculate the far-field approximation of ${\bar{G}}^{\mu}_{\tilde{\theta}\;\nu}$. Starting from the second expression in Eqs. (43) and choosing the coordinate system ${\cal S}$ previously introduced in Fig. 7, ${\bar{G}}^{\mu}_{\tilde{\theta}\;\nu}$ can be rewritten as

[TABLE]

where we define the following integrals:

[TABLE]

which we subsequently separate according to $\alpha=0$ and $\alpha=k=1,2$. Rewriting the angular integration in terms of Bessel functions and employing the expression (89) for the Sommerfeld integral with $Z=|z|+|z^{\prime}|$, we obtain

[TABLE]

Next, we implement the far-field conditions (86), yielding

[TABLE]

in the limit $\rho\rightarrow\infty$ and $z\rightarrow\infty$. Then, the following product reduces to

[TABLE]

in each of the integrands of Eqs. (92) and (93). This factor is a rapidly oscillating function of $k_{\perp}$ that allows one to apply the stationary phase approximation to evaluate the dominant contribution Chew ; Chew1 . Recalling that $k_{z}=\sqrt{\tilde{k}_{0}^{2}-k_{\perp}^{2}}$, the stationary phase condition is

[TABLE]

with $r=\sqrt{\rho^{2}+z^{2}}$. Moreover, in Eqs. (92) and (93) we can estimate $\tilde{k}_{0}^{2}-k_{\perp}^{2}=k_{z}^{2}\simeq(k_{z})_{s}k_{z}$ around the point where the phase is stationary. Using again the Sommerfeld identity (89) with ${\cal R}={\tilde{R}}$ and $Z=|z|+|z^{\prime}|$, we arrive at

[TABLE]

where

[TABLE]

Finally, we complete the far-field approximation in Eqs. (97) and (98) by writing

[TABLE]

with $\mathbf{n}_{\perp}=\sin\theta(\cos\phi,\,\sin\phi,\,0)$. The results are

[TABLE]

where we have dropped terms of $\mathcal{O}(r^{-2})$ and higher.

Next we comment on a technical point, which fortunately is not relevant to our purposes. On one hand, the integrand in Eq. (92) has poles at $k_{\perp}=\pm{\tilde{k}}_{0}$. At the same time, the stationary point obtained in Eq. (96) indicates that $(k_{\perp})_{s}={\tilde{k}}_{0}\sin\theta$. Both values coincide for $\theta=\pm\pi/2$ indicating the factor of the exponential is not a smooth function around the stationary point now. In this way, the stationary phase calculation requires some modifications, which we do not pursue in this work. The reason we do not require such improvements is because the VC cone condition demands $0<{\theta}_{\rm C}<\cos^{-1}(1/n)$, which is always far away from the dangerous point $\theta=\pi/2$, unless $n$ is extremely large. In other words, our calculation in the stationary phase method is still a good approximation for the far-field behavior of our GFs when $\theta<\pi/2$. An analogous situation occurs when dealing with the steepest descent method WAIT ; VANDER ; CLEMM .

At this point and as a matter of consistency, we show that the stationary phase method and the steepest descent method lead to the same results, under a similar restriction for ${\theta}_{\rm C}$. Let us recall that the latter method Chew ; Mandel allows the calculation of the leading contribution to the integral

[TABLE]

when a rapidly varying exponential factor multiplies the function $f(t)$, which must have a smooth behavior in the region close to the stationary phase point, yielding the result

[TABLE]

where $t_{0}$ is such that $h^{\prime}(t_{0})=0$. From Eqs. (92) and (93), we extract the required integral as

[TABLE]

whose asymptotic behavior for $\tilde{k}_{0}^{2}>k_{\perp}^{2}$ is

[TABLE]

Since the rapidly oscillating phase in Eq. (106) coincides with that of Eq. (95) in the far-field regime, the stationary point of $K_{1}(\mathbf{x};\omega)$ is also given by the conditions in Eq. (96), yielding

[TABLE]

as the phase at the stationary point. After this, we require computation of the second derivative of the phase at the stationary point in the far-field regime, which is

[TABLE]

Substituting Eqs. (107) and (108) in Eq. (104), we obtain that the leading-order term calculated by this method is

[TABLE]

where we have disregarded ${\pi}/{4}$ in front of $\tilde{k}_{0}r\rightarrow\infty$ and we have replaced ${\tilde{R}}$ by its expression (100). Inserting Eq. (109) in the expressions (92) and (93), we obtain the results given in Eqs. (101) and (102) by the stationary phase method. Therefore, putting together Eqs. (101) and (102) in Eq. (90), we find the final expression for ${\bar{G}}^{\mu}_{\tilde{\theta}\;\nu}(\mathbf{x},\mathbf{x}^{\prime};\omega)$ as

[TABLE]

where ${s}_{\alpha}=(1/n,\mathbf{\hat{n}})$.

The calculation of $\bar{G}_{\tilde{\theta}^{2}\;\nu}^{\mu}$ proceeds along similar lines and the final result is summarized in the equations

[TABLE]

where

[TABLE]

We have also verified that the steepest descent method Chew ; Mandel gives the same results as the stationary phase approximation in the calculation of the required rapidly oscillating integrals for this contribution.