Radiation of a short linear antenna above a topologically insulating half-space

The topological magnetoelectric effect (TME) is a unique macroscopic manifestation of quantum states of matter possessing topological order and it is described by axion electrodynamics. In three-dimensional topological insulators, for instance, the axion coupling is of the order of the fine structure constant, and hence a perturbative analysis of the field equations is plenty justified. In this paper we use Green's function techniques to obtain time-dependent solutions to the axion field equations in the presence of a planar domain-wall separating two media with different topological order. We apply our results to investigate the radiation of a short linear antenna near the domain-wall.


I. INTRODUCTION
Topological materials are not distinguished by local ordered parameters, but by topological quantum numbers which could be described and computed using topological field theories in the low energy limit.The best studied of these are the topological insulators (TIs), where the corresponding bulk band gap is similar to an ordinary insulator, but the boundary of the material hosts gapless states which are protected by the topological numbers of the bulk geometry [1].The experimental discovery of TIs in HgTe quantum wells [2] and stoichiometric crystals [3][4][5] encourages further investigations regarding its exotic properties.
Apart from their spectroscopic distinguishing features, TIs have unconventional electromagnetic response described by axion electrodynamics, which is defined by the standard Maxwell Lagrangian plus the additional θ-term (in SI units) L θ = (α/π) ϵ 0 /µ 0 θ E • B, where E and B are the electromagnetic fields, α is the fine structure constant, and θ is an angular variable known in particle physics as the axion field or axion angle [6].The axion field was originally proposed to solve the strong-CP problem in QCD [7][8][9], but its weak interaction strength with standard model particles and an elegant production mechanism in the early universe make it a promising candidate to explain dark matter as well [10][11][12].In the context of TIs, time-reversal (TR) symmetry indicates that θ = 0, π (mod 2π), and hence the θ-term has no effect on Maxwell equations in the bulk.The nontrivial topological property, a surface half-integer quantum Hall effect, manifests only at the boundary between two manifolds with different topologies, in which case θ is not globally constant.Therefore, from the macroscopic point of view, the axion angle θ can be viewed as an effective parameter characterizing the media, in the same footing as its permittivity ϵ and permeability µ.
Time-independent solutions to the axion field equations have been elaborated by using the electric-magnetic duality of Z 2 TIs [15] as well as Green's function techniques [17,19].However, general solutions to the time-dependent case are still under investigation [40][41][42].In this work we aim to fill in this gap by computing, in a perturbative fashion, solutions to the timedependent axion field equations in the presence of a planar domain wall separating two media with different topologies.The validity of a perturbative scheme for the solution of the field equations is as follows.In the high-energy context the validity is plenty justified: the axion-to-two-photon coupling g θγγ which converts photons into axions is expected to be very small [43].In the condensed matter framework, the smallness of the θ-term is controlled by the coupling αθ, which is small since α ≈ 1/137 and θ = π for a TI with a single Dirac fermion on its surface.In the hypothetical situation of a TI with a large number of Dirac fermions on its surface, our perturbative approach ceases to be valid as αθ approaches unity.In view that exact solutions to the static axion field equations have been reported in the literature, we shall first validate our perturbative scheme by reproducing well-known results in this case.Next we will safely extend our analysis to the time-dependent case.Finally, we apply our results to study the radiation of a short linear antenna near the domain-wall.
The paper is organized as follows.In Sec.II we briefly review the basics of axion electrodynamics.In Sec.III we solve, in a perturbative fashion, the static field equations and verify its consistency with previously reported results for the Green's function for two ponderable topological insulators separated by a planar interface.Section IV is devoted to the evaluation of the time-dependent Green's function by assuming equal optical properties in the two media.This configuration allow us to isolate the topological effects.The details of technical calculations are relegated to the appendices.Finally, in Sec.V we apply our result to investigate the radiation of a short linear antenna near a θ domain-wall.

II. TOPOLOGICAL FIELD THEORY
The effective field theory governing the electromagnetic response of topological insulators, independently of the microscopic details, is defined by the action (in SI units) [44]: where α = e 2 /(4πϵ 0 ℏc) ≃ 1/137 is the fine structure constant, ϵ and µ are the permittivity and permeability of the material, respectively, and θ is the topological magnetoelectric polarizability (axion field).Under periodic boundary conditions, the partition function and all physical quantities are invariant under shifts of θ by any multiple of 2π.Since E • B is odd under TR, the only symmetry allowed values of θ are 0 or π (modulo 2π).This leads to the Z 2 classification of three-dimensional TR invariant TIs.Clearly, the θ-term in Eq. ( 1) has no effect on Maxwell equations in the bulk, since it is a total derivative term.The nontrivial topological property, a surface half-integer quantum Hall effect, manifests only when a TR-breaking perturbation is induced on the surface to gap the surface states, thereby converting the material into a full insulator.This can be achieved by introducing magnetic dopants to the surface [45] or by the application of an external static magnetic field [46].In this situation, θ is quantized in odd integer values of π, where the magnitude and sign of the multiple is controlled by the strength and direction of the TR-breaking perturbation.
The electromagnetic response of a system in the presence of the axionic term is still described by the ordinary Maxwell equations, but the constituent relations which define the electric displacement D and the magnetic field H acquire an extra term proportional to θ, i.e.D = ϵE + α(θ/π) ϵ 0 /µ 0 B and H = µ −1 B − α(θ/π) ϵ 0 /µ 0 E. Here, the θ-dependent terms, which arise from the axionic term in Eq. ( 1), describe the topological magnetoelectric effect.The frequency components of the electromagnetic fields (introduced as usual: by making a Fourier decomposition in time of the sources and fields) obey the inhomogeneous equations while the homogeneous equations are untouched due to gauge invariance of the action (1): ∇ × E(r, ω) = iωB(r, ω).
One can further verify that the divergence equations ( 2) and ( 4) are subsumed within the curl equations ( 3) and ( 5), which can be verified by taking the divergence of the latters.We observe that equations ( 2)-( 3) depend only on the space gradient of the axion field.So, the θ-term does not modify the field equations in the bulk of topological insulators.However, for a TI in contact with a conventional insulator, as shown in Fig. 1, the axion field is piecewise constant across the interface Σ and then ∇θ(r) = (θ 1 − θ 2 )δ(Σ)n, where n is the outward unit normal to Σ.These are the equations we aim to solve in a perturbative fashion along this paper.Our general approach will be to first Taylor expand the electromagnetic potentials in powers of the fine structure constant and then solve recursively the resulting equations.Before embarking to the solution of the time-dependent field equations, we first tackle the time-independent case and demonstrate consistency of our results with previously reported ones.

III. SOLUTION TO THE STATIC AXION FIELD EQUATIONS
Let us consider the geometry shown in Fig. 1 with the interface Σ being the plane z = 0.The half-space z < 0 is occupied by a topological medium with dielectric constant ϵ 2 , magnetic permeability µ 2 and axion angle θ 2 , whereas the half-space z > 0 is occupied by a medium with dielectric constant ϵ 1 , magnetic permeability µ 1 and axion angle θ 1 .Then the inhomogeneity in the permittivity, permeability and the axion field are limited to a finite discontinuity across the planar surface z = 0. To proceed further we write the field equations ( 2) and (3) in terms of the electrostatic potential ϕ(r) and vector potential A(r).In the Coulomb gauge, the static electromagnetic potentials satisfy where α = α(θ 1 −θ 2 )/π.To obtain a perturbative solution for the electromagnetic potentials in the presence of arbitrary external sources, we introduce the expansions ϕ(r , where ϕ (n) and A (n) (r) are understood to vanish as αn .We now substitute the power expansion into equations ( 6) and (7), and comparing coefficients of each power of α results in an infinite series of simultaneous equations.The zeroth-order equations are simply the Maxwell equations in the absence of the θ-term, The n-order equations, which are the main concern in this work, are The process of solving these equations begins with the solution to the zeroth-order equations (8), which involves the wellknown electro-and magnetostatic Green's functions.We now adopt the following notation.Since translational invariance is broken in the z direction due to the presence of the interface, it is natural decompose any vector C into components parallel and perpendicular to the normal direction êz , i.e.C z = C z êz and C ⊥ = C x êx + C y êy , respectively.So, the general solution for the electrostatic potential is where the electric Green's function is given by in which r ⊥ = xê x + yê y and r ′ ⊥ = x ′ êx + y ′ êy [47].The solution to the magnetostatic problem is a bit more subtle.Decomposing the current density and the vector potential into their components parallel and perpendicular to êz , i.e.J = J ⊥ +J z and A = A ⊥ + A z , respectively, the general solution for the vector potential, in the half-space z > 0, can be written as while for the region z < 0 the solution is simply The corresponding magnetic Green's function is given by [47] where χ = +1 (−1) for ξ = ⊥ (z).With all of the above we are ready to solve equations (9).This is a simple but not straightforward task, so we relegate the details to Appendix A. The final results, up to second order in the parameter α, are where the indices i, j range over the transverse coordinates x, y, and with the definition Clearly, the electric (vector) potential is even (odd) under time-reversal operation, as expected.

Image magnetic monopole effect
As a consistency check of our results, we shall verify the image magnetic monopole effect, which has been extensively discussed in the literature by using different approaches, for example, the method of images [14], the electric-magnetic duality of TIs [15] and Green's function techniques [19,20].This elusive effect consists in the appearance of an image magnetic monopole when an electric charge is placed at a certain distance from the TI's surface.
We now consider for the charge density a pointlike electric charge q at z = z 0 .Equations ( 14) and ( 15), together with the Green's functions ( 16)-( 18), give the electric and vector potentials.For z > 0, the electric potential can be interpreted as due to two pointlike electric charges, one of strength q at z = z 0 , and the other, the image charge, of strength at −z 0 .For z < 0 only one pointlike electric charge appears, of strength q ′ +q located at z 0 .Unremarkably, the above expression reduce to the standard textbook one in the nontopological limit (α = 0) [47].Due to the magnetoelectricy of TIs, a magnetic field is induced.The corresponding vector potential is found to be which we identify as the vector potential of a straight vortex line or Dirac string over the z axis.For z > 0 (z < 0) the string is over the negative (positive) z axis ending at a magnetic monopole at −z 0 (z 0 ).In other words, the magnetic field for z > 0 can be interpreted as that of a magnetic monopole of strength at the image point −z 0 .Similarly, for z < 0 the magnetic field can be understood as originated by a monopole of strength −g located at z 0 .Our results are in a agreement with the ones reported in Refs.[14,15,19,20].This can be further verified by power expanding the image electric and magnetic charges up to order α2 , thus yielding our results ( 20) and ( 22), respectively.

IV. TIME-DEPENDENT SOLUTIONS TO THE AXION FIELD EQUATIONS
In this Section we consider the same configuration as in Sec.III, but we shall investigate now time-dependent solutions to the field equations ( 2)- (5).In order to isolate the contribution arising from the surface Hall effect, we restrict ourselves to the case in which the permittivities are the same, i.e. ϵ 1 = ϵ 2 ≡ ϵ, but the topological magnetoelectric polarizability is a piecewise constant function: θ 1 for z > 0 and θ 2 for z < 0. We assume the sources located in the region z > 0.
while the n-order potentials are determined by The first step towards solving the higher-order equations ( 27)-( 28) is to solve the zeroth-order equations ( 25)- (26).This requires the well-known frequency-dependent Green's function of electromagnetism [47]: such that With these results we can evaluate now equations ( 27) and (28).In this case we solve the system up to first order in the parameter α.As shown in the Appendix B 1, for the vector potential one gets: and for the scalar potential one spot the identity i ω v 2 ϕ(r, ω) = ∇ • A(r, ω).Here, where A detailed derivation of this expression can be found in the Appendix B 2. It is worth to mention that the exact electromagnetic potentials (and hence the fields) must be calculated using Eq. ( 33), which unfortunately cannot be evaluated in an analytical FIG. 2. Short linear antenna (of length 2d and frequency ω) at a distance z0 from the planar interface between two topologically insulating media.
fashion.However, for many practical calculations, one can obtain an approximated expression valid for In the Appendix B 3 we derive, by two different methods, the aforementioned expression.The approximation relays in the fact that the integrand in Eq. ( 33) is a rapidly oscillating function for large |r ⊥ − r ′ ⊥ | and |z| + |z ′ |, and hence the leading contribution can be obtained by using, for example, the stationary phase approximation.The final result is Once we have found the indexed Green's function, we are ready to tackle problems involving different time-dependent sources near the domain wall.

V. RADIATION OF A SHORT LINEAR ANTENNA NEAR A θ DOMAIN-WALL
The problem of a short linear antenna near a dielectric half-space has attracted great attention for over the years due to its importance to many practical applications, in particular, in the context of THz applications [48], near-field optics [49], plasmonics [50] and nanophotonics [51].The exact solution to this problem was given first by Sommerfeld in 1909 [52].Since then numerous investigations have been carried out in this subject.The main goal of this Section is to address the problem of the radiation pattern produced by a short linear antenna placed near the boundary between two θ-half spaces.To this end, we shall employ the theory developed in the previous Section.
The geometry of the problem is shown in Fig. 2. The left half-space z < 0 is occupied by a topological insulator with axion angle θ 1 , while the right half-space z > 0 is occupied by a TI with axion angle θ 2 .We consider a simple model for a transmitting antenna oriented parallel to the interface, say along the the x-axis and placed at a distance z 0 from the interface, with uniform current flowing along a conductor of length 2d.We assume that this current varies sinusoidally with time as where I 0 and ω 0 denotes the current amplitude and angular frequency, respectively, the δ's are Dirac deltas and Θ is the Heaviside step function.Because the current is discontinuous at the ends, charge must be deposited there, being of opposite sign at each end, thus forming an electric dipole.The exact zeroth-order scalar and vector potentials can be calculated in a simple fashion and can be found in a variety of textbooks on electrodynamics [47].The only nonzero component of the vector potential (which is in the same direction as the current) in the time-domain is where r ′ (x ′ ) = x ′ êx + z 0 êz .The result for the scalar potential in the time-domain is where r ± = r ′ (±d).One can easily verify that these potentials satisfy the Lorentz gauge condition.In order to evince the impact of the axion field, we first evaluate numerically the electromagnetic potentials and fields at first order in θ by using Eq. ( 31) together with the exact Green's function (32) with G ω given by Eq. (33).To this end we introduce the dimensionless quantities X = x/z 0 , Y = y/z 0 and Z = z/z 0 .Figure 3  The former exhibits the usual current-encircling behavior dictated by the Ampère law, while the latter exhibits the expected results, whose filed lines invert their direction periodically.In order to show the full electromagnetic fields, including the effects of the magnetoelectric coupling, we have to introduce further dimensionless quantities besides the coordinates (X, Y, Z), namely, T = t/τ , D = d/z 0 , Λ = ω 0 z 0 /c and τ = z 0 /c.This is so because the appearance of an additional length scale in the system: the antenna-surface distance z 0 .In Fig. 4 we show the total electromagnetic fields for θ = π, which is a typical value for topological insulators.We observe that the cross-sections X = 0 exhibit the same behavior as in the case θ = 0, for both, the electric (panel (a)) and magnetic (panel (c)) fields; however, as we move from X = 0 the field lines deviates from the nontopological case.Panels (b) and (d) show the Y = 0 cross-section, in which case we observe a considerably different behavior, which cannot be interpreted in a simple manner.Indeed, the charge and current densities sourcing these fields are far from trivial, but they can be understood in terms of the half-integer Hall current at the interface.The normalized Hall current j Hall (in units of σ Hall E 0 , with σ Hall = e 2 /2h) is where with Now we are interested in the fields that survive at large distances from the source, i.e. in the radiation zone, which requires approximating the potentials for r ≫ λ and r ≫ d, where λ = 2πc/ω is the wavelength.To this end we use the expansions , where n = r/r.All in all, the zeroth order potentials (up to order 1/r) become Using these potentials we find the electromagnetic fields in the radiation zone: Following the same program one can derive expressions for the electromagnetic potentials at first-order.The details are relegated to the Appendix C. Using Eq. ( 31) one finds that the nonzero components of the vector potential in the time-domain at first-order become where The first-order scalar potential in the time-domain is In the radiation zone these potentials can be further simplified.Indeed, one realizes than only the component A y of the vector potential contributes to radiation (to be precise, the term proportional to k 2 0 in Eq. ( 43)), while the others decay faster than 1/r.

Using the expansion |r
with r * = r, we get (r, t) = v y r A (1)  y (r, t).(45) Keeping terms of the order 1/r the electromagnetic fields become As shown in Ref. [17], energy conservation is unchanged in the presence of the θ term, i.e. given the Poynting vector S and energy density u defined as usual S = E × H and u = 1 2 (E • D + B • H), respectively, one realizes that ∂ t u + ∇ • S = 0. Momentum conservation is modified by a surface term due to the self-induced charge and current densities arising there, but it does not concern us in this paper.The angular distribution of the radiated power takes the usual form dP dΩ = v µ r 2 |B| 2 , with however a total magnetic field which acquires a contribution arising from the θ-term.Keeping the lowest order in the magnetoelectric polarizability α and averaging (in time) over a complete cycle, the angular distribution of the radiated power by the linear antenna becomes The first term corresponds to the usual radiated power by a linear antenna embedded in an infinite material media described by the permittivity ϵ and permeability µ.The second term, which is proportional to the surface Hall conductivity, exhibits an interesting divergent behavior near the interface, as is depicted in Fig. 5.The physical origin of this term is understood from the surface Hall current and charge densities at the interface, which generates electromagnetic field which contribution to radiations in the near zone to the interface.

VI. SUMMARY AND DISCUSSION
Topological phases are an emerging class of materials which have attracted much attention in condensed matter physics.The remarkable microscopic properties encoded in the band structure of topological phases produces exceptional macroscopic electromagnetic response, which is described by a topological field theory akin to axion electrodynamics.In general, a linear magnetoelectric media can be described by adding the term ∼ θ ij E i B j to the free enthalpy of the system, where θ ij is the magnetoelectric tensor.In this frame, axion electrodynamics is a particular case when the ME tensor is of the form θ ij = θδ ij .
In this paper we derived, in a perturbative fashion, solutions to the time-dependent axion electrodynamics field equations for a picewise constant axion angle.The validity of the perturbative program relies in the fact that the axion coupling, in the case of topological insulators, is proportional to the fine structure constant.In the particle physics context this condition is vastly fulfilled, since it is expected to be very small.The difference is that in the former case the presence of boundaries is a natural task, while in the latter the axion is assumed to be distributed over the whole space.
To validate our perturbative scheme, we first tackle the static axion field equations and derive the components of the Green's matrix up to second order in the axion coupling.We show that our results are consistent with the previously reported exact Green's matrix when Taylor expanded up to second order in the coupling.We discuss in passing the image magnetic monopole effect, which is a salient manifestation of the unconventional magnetoelectric response of topological insulators.In the timedependent case, the axion field equations become more involved and no analytical results have been reported yet, although there are some first attempts towards this direction.we aim to fill in this gap by applying the same perturbative program to the axion field equations.For the sake of simplicity, in this paper we choose that the bulks have the same optical properties, namely, ϵ 1 = ϵ 2 ≡ ϵ and µ 1 = µ 2 ≡ µ, which is made to highlight and isolate the effects of the axion coupling.In a forthcoming paper we will relax this condition and generalize our results to the case of two topologically insulating media characterized by different nontrivial optical parameters.This generalization is of course far from trivial.To illustrate the importance of the axion angle in the time-dependent case considered here, we investigate the radiation of a short linear antenna above a topologically insulating domain-wall with different values of θ.
We first computed numerically the electromagnetic fields due to the linear antenna by using the exact Green's function and verify that they exhibit the expected behavior.We interpret the induced fields in terms of the Hall current at the interface between the two half-spaces characterized by different θ-angles.Next, to investigate the effects of the θ-term in the power radiated, we use the asymptotic expression for the Green's function and evaluate the fields in the radiation zone.These field allow us to compute the time-averaged angular distribution of the radiated power, which as expected, is dominated by the usual Maxwellian term, but becomes corrected due to the θ-term.Such correction is quite interesting, since exhibits a divergent-type behavior near the interface, which arises from the field produced by the surface charge and current densities there.This is a distinctive behavior with respect to the usual radiation pattern.
In summary, a general solution to the static axion field equations has been previously elaborated in terms of Green's function techniques; however, the time-dependent case has received less attention mainly due to its complexity.In this paper we introduced a perturbative scheme which was successfully validated with the static theory and then we extended the analysis to the time-dependent case.We obtained general expressions for the corresponding matrix of Green's function which we also approximated in the far-field regime.We expect our method and results will be of considerable relevance and that they may constitute the basis for numerous other researches which include more intricate model of antennas and the radiation of particles near the θ-interface.
Using the electric (11) and magnetic ( 13) Green functions we realize that To carry out integration, we use the following integral representation: where Σ k ⊥ is the two-dimensional space k x -k y .Using this representation twice in Eq. (A6) we get Integration with respect to the coordinates r ′′ ⊥ produces a Dirac delta on the momenta k ⊥ and k ′ ⊥ , i.e.
where we used the properties of the Dirac delta to perform integration over k ′ ⊥ .The resulting double integral becomes easier to perform if we express the area element in polar coordinates, i.e. d 2 k ⊥ = k ⊥ dk ⊥ dφ, and choose the k x axis in the direction of the vector The parentheses in this equation enclose an integral representation of the Bessel function J 0 .The resulting integral is and taking the primed derivative This integral is well known [53].The final result is then 2. Second-order Green's functions We now turn to the second-order corrections to the GFs.Integration of Eq. ( 9) produces where Using the first-order Green's function given by Eq. (A13), we verify that ∂ ′′ j G (1) i0 (r ′′ , r ′ ) for j = 1, 2, and hence the double primed derivatives can be get out of the integrals as a minus primed derivative, as before.Performing the integration with respect to z ′′ we get Using the integral representations (A4) for the 0i-components and (A5) we find and substituting the electric (11) and magnetic (13) Green functions we get where we have defined the symmetric matrix I = (I ij ) ∈ R 2×2 , with components To evaluate this expression, we use thrice the integral representation (A7), i.e.
Integrating with respect to r ′′ ⊥ and r ′′′ ⊥ produces two Dirac deltas in momentum: where we have used the properties of the Dirac delta.We now introduce the polar coordinates (k x , k y ) = k ⊥ (cos φ, sin φ) and choose the vector r ⊥ − r ′ ⊥ along the k x direction.So, the angular integral yields This integral does not converge as is, but its derivatives does.Taking one of the derivatives we find To perform this integral in a convergent manner, we use the following trick: ∞ Λ e −k ⊥ λ dλ = e −k ⊥ Λ /k ⊥ , and therefore and this integral is well-known; see, for example, Ref. [53].The result is where we have defined We now perform son manipulations in Eqs.(B1) and (B2).Taking the divergence of Eq. (B2) we obtain where we used that ϵ zij ∂ i ∂ j = 0. So, one we have computed the vector potential, the scalar potential becomes fully determined.
In order to manipulate Eq. (B2) we use the local charge conservation: ∇ • J = iωρ.Substituting the charge density in Eq. (B2) produces Integrating by parts and using the identity which can be written as where This expression corresponds to Eq. ( 32) of the main text.
In the following we derive an exact expression for the function G ω (r, r ′ ), which cannot be expressed in term of known functions.We also obtain an analytical expression in the radiation zone (far field approximation) by using two different methods.The integral we have to evaluate is given by Eq. (B3).Substituting the frequency-dependent Green's function (29) the integral becomes: To evaluate this integral we now use the Sommerfeld identity: where Σ k ⊥ is the two-dimensional space k x -k y and r = r ⊥ + |z|ê z .Using this representation twice in Eq. (B9) we get Integration over r ′′ ⊥ produces a Dirac delta in momentum, which allow us to evaluate one of the momentum integrals.We obtain: As before, we now introduce the polar coordinates (k x , k y ) = k ⊥ (cos φ, sin φ) and choose the vector r ⊥ − r ′ ⊥ along the k x direction.So, the angular integral yields This is an exact expression for the required function.It corresponds to Eq. ( 33) of the main text.

Far-field approximation
Here we shall approximate the integral (B13) in the radiation zone, i.e. for ρ, |z| ≫ 1.For crosscheking our result, we compute the integral by two different methods.One of them used the Sommerfeld integral (Method 1) and the other is the stationary phase approximation (Method 2).Method 1.We start by using the identity q −1 = ∞ 0 dξ e −qξ in Eq. (B13), for q Using the Sommerfeld integral (B10) this expression can be written in the simple form: where Next we shall approximate this integral in the far-field region.To this end we use the fact that R where ζ = k 0 R 0 and z = (|z| + |z ′ |)/R 0 .Numerical evaluation of this integral for ζ ≫ 1 and z ≪ 1 reveals that only small values of χ contribute to the integral.So we can expand the root and retain solely the leading order, i.e.
This result corresponds to Eq. ( 34) of the main text.
Method 2. Since the Bessel function and the complex exponential in the integrand of Eq. (B13) are rapidly oscillating functions in the far-field region, one can use the stationary phase approximation to evaluate the integral.To this end we first use the identity J 0 0 (x) + H (2) 0 (x)] together with the reflection formula H 0 (e iπ x) = −H (2) 0 (x), where H 0 and H (2) 0 are the Hankel functions.In this case, the integral (B13) can be written as We now implement the far-field approximation into the Hankel function.For x ≫ 1 we have H We observe that, in the far-field, the exponential factor is a rapidly oscillating function of k ⊥ and hence we can apply the stationary phase approximation to evaluate the dominant contribution.A simple calculation reveals that the stationary phase point locates at Finally, we use the well-known asymptotic formula R g(x)e iλf (x) dx ≈ g(x 0 )e iλf (x0)+i π 4 sgn(f ′′ (x0)) 2π λ|f ′′ (x 0 )| , (B22) where λ ≫ 1 and f ′ (x 0 ) = 0, to obtain which is again the result of Eq. (B18).This validates our expression for the far-field approximation.
Appendix C: Evaluation of the first-order potentials for a short linear antenna Our starting point here is Eq. ( 31) for the components of the vector potential.The first order contribution in the frequencydomain is read off as where G ω (r, r ′ ) is the corresponding Green's function.The potential in the time-domain is computed by inverse Fourier transforming the above expression.We obtain: The current source (35)  where r ′′ = x ′ êx + z 0 êz .Since we are interested in the far-field domain, we use the Green's function (34).In such case the vector potential becomes where r ′ i (x ′ ) = x ′ êx − z 0 êz is the image point of r ′ (x ′ ) and r * = xê x + yê y + |z|ê z .Equations ( 42) and ( 43) of the main text are derived from this result by taking i = x and i = y, respectively.Now, to derive the scalar potential, we appeal to the Lorentz gauge condition.Taking the divergence of the vector potential we obtain wherefrom we find which, after taking the y-derivative, corresponds to Eq. ( 44) of the main text.

FIG. 1 .
FIG.1.Geometry of two semi-infinite topologically insulating media separated by the surface Σ.The media have different magnetoelectric polarizabilities, permittivities and permeabilities, as depicted.

FIG. 3 .FIG. 4 .
FIG. 3. Illustration of cross-sections of the normalized electric field [panels (a) and (b)] and magnetic field [panels (c) and (d)] generated by a short linear antenna in vacuum at z0.
shows the electric and magnetic fields in the case θ = 0, i.e. in the absence of the topological magnetoelectricity.Panels (a) and (b) show the electric field lines (in units of E 0 = µIv 4πz0 ) on the planes Y = 0 and X = d/Z 0 (where the end-points of the antenna are located) for ωt = π/4, respectively.They exhibit the expected behavior: dipole-like in panel (a) and almost point-like in panel (b).Panels (c) and (d) show the magnetic field lines (in units of B 0 = E 0 /v) on the planes X = d/z 0 and Y = 0 for ωt = π/4, respectively.

FIG. 5 .
FIG. 5. Angular distribution of the radiated power by a thin linear antenna near a θ domain-wall.The position of the interface is evident.