Absence of induced magnetic monopoles in Maxwellian magnetoelectrics

The electromagnetic response of topological insulators is governed by axion electrodynamics, which features a topological magnetoelectric term in the Maxwell equations. As a consequence magnetic fields become the source of electric fields and vice-versa, a phenomenon that is general for any material exhibiting a linear magnetoelectric effect. Axion electrodynamics has been associated with the possibility to create magnetic monopoles, in particular by an electrical charge that is screened above the surface of a magnetoelectric material. Here we explicitly solve for the electromagnetic fields in this geometry and show that while vortex-like magnetic screening fields are generated by the electrical charge their divergence is identically zero at every point in space which implies an absence of induced magnetic monopoles. Nevertheless magnetic image charges can be made explicit in the problem and even if no bound state with electric charges yielding a dyon arises, a dyon-like angular momentum follows from our analysis. Because of its dependence on the dielectric constant this angular momentum is not quantized, which is consistent with a general argument that precludes magnetic monopoles to be generated in Maxwell magnetoelectrics. We also solve for topologically protected zero modes in the Dirac equation induced by the point charge. Since the induced topological defect on the TI surface carries an electric charge as a result of the axion term, these zero modes are not self-conjugated.

The electromagnetic response of topological insulators is governed by axion electrodynamics, which features a topological magnetoelectric term in the Maxwell equations. As a consequence magnetic fields become the source of electric fields and vice-versa, a phenomenon that is general for any material exhibiting a linear magnetoelectric effect. Axion electrodynamics has been associated with the possibility to create magnetic monopoles, in particular by an electrical charge that is screened above the surface of a magnetoelectric material. Here we explicitly solve for the electromagnetic fields in this geometry and show that while vortex-like magnetic screening fields are generated by the electrical charge their divergence is identically zero at every point in space which implies an absence of induced magnetic monopoles. Nevertheless magnetic image charges can be made explicit in the problem and even if no bound state with electric charges yielding a dyon arises, a dyon-like angular momentum follows from our analysis. Because of its dependence on the dielectric constant this angular momentum is not quantized, which is consistent with a general argument that precludes magnetic monopoles to be generated in Maxwell magnetoelectrics. We also solve for topologically protected zero modes in the Dirac equation induced by the point charge. Since the induced topological defect on the TI surface carries an electric charge as a result of the axion term, these zero modes are not self-conjugated.

I. INTRODUCTION
A remarkable feature of three-dimensional topological insulators (TIs) is their so-called magnetoelectric (ME) effect, a collection of phenomena where magnetic fields become the source of electric fields and vice-versa [1]. This topological electromagnetic response is governed by so-called axion electrodynamics, which features a magnetoelectric term L a = αθ/(4π 2 )E · B in the Lagrangian density L a , with electric and magnetic fields E and B respectively, θ a 2π-periodic parameter and α the finestructure constant. In a topological insulator θ is a parameter that follows from the band structure topology, being given by a Berry non-Abelian flux in the Brillouin zone [1,2]. By symmetry the magnetoelectric coupling term is actually present in any material that exhibits a linear magnetoelectric effect -induction of magnetization by an electric field or of electric polarization by a magnetic field. However, in ordinary magnetoelectric materials such as Cr 2 O 3 , BiFeO 3 , and GdAlO 3 the magnetoelectric coupling constants are quite small [3]. The topological ME effect has been recently measured using Faraday and Kerr rotation [4,5], which were shown to be quantized according to the prediction of axion electrodynamics of TIs.
A number of further interesting consequences of the axion term have been predicted, for instance that a cylindrical TI becomes electrically polarized under an applied magnetic field parallel to the cylinder symmetry axis [6]. An interesting possible experimental setup exploring this effect is a flux tube piercing the interior of a TI, perpendicular to its surfaces [7]. If the surfaces are coated with thin film ferromagnets with opposite magnetizations, the surface states become gapped and a topological electro-magnetic response ensuing the axion term in the Lagrangian occurs. In this scenario the cylinder becomes an Aharonov-Bohm flux tube and an electrical polarization is induced leading to fractional charges ±e/2 on the top and bottom surfaces, respectively [7]. In the case of a magnetic vortex that enters from a superconductor (SC) into a time-reversal invariant TI, it was shown that the vortex induces a charge of e/4 [8,9] at the SC-TI interface. In this situation also the vortex angular momentum, which determines the vortex statistics, is fractional [9]. The emergence of fractional charges is reminiscent of the Witten effect [10], which predicts that the axion term causes electric charge fractionalization in the presence of magnetic monopoles. When real magnetic monopoles were to be present in an axion magnetoelectric, fractional electric charges would occur not only at surfaces but also in the bulk of a magnetoelectric, since in this case the Maxwell equations are modified despite the axion Lagrangian L a being a total derivative [11].
In this context it is highly interesting that the presence of a magnetoelectric term in the Maxwell Lagrangian has been associated with the possibility of creating magnetic monopoles. In particular, the situation has been considered in which a magnetic monopole emerges from the screening of an electrical charge that is situated outside a magnetoelectric material, at a certain distance d from its the surface [12,13], see Fig. 1. This would be quite remarkable as it would imply that the condition that the magnetic field be divergence-free is lifted by the axion term in one way or another. Here we revisit this very well defined geometry and determine the unique solution for the electric and magnetic fields by direct evaluation, without resorting to an image charge construction [12] or Green function formalism [14] and determine also the angular momentum of this dyon-like object carrying both electric and magnetic charge. The resulting divergence of the magnetic field vanished at every point in space. Induced magnetic monopoles are thus absent and instead the electrical charge generates a magnetic vortex structure near the magnetoelectric surface which in turn generates magnetic screening fields in all of space. In the limit that the electrical charge is placed at the TI surface (d → 0) it almost behaves like a magnetic monopole, but actually corresponds to a point vortex, or Pearl vortex [15], still satisfying the local constraint ∇ · B = 0. Although the solution can be cast in terms of image electric and magnetic charges, the latter are not related by a Dirac duality quantization characteristic of magnetic monopoles and hence cannot be interpreted as such. We will show that despite this difficulty dyon-like behavior occurs. Dyons are dipoles constituted of an electric and a magnetic charge and are well known to have an angular momentum [16] with the following properties: (i) its value is independent of the separation between the electric and magnetic charges, (II) it is nonzero even if the dyon is at rest, and (iii) it is quantized of the Dirac duality relation eg/c = n /2, n ∈ Z. For our problem of a point charge a distance d apart from the TI surface, the resulting angular momentum fulfils properties (i) and (ii), but not (iii), which is signalled by the angular momentum explicitly depending on the dielectric constant of the TI. Insisting the angular momentum be quantized, would imply the (static) dielectric constant be negative, violating the inequality > 1. This dependence on prevents the interpretation of this dyon-like object as an anyon, a quasi-particle having fractional statistics as quantum statistics is a universal property of (quasi-)particles and should not depend on the details of the medium they are embedded in.
The outline of the paper is as follows. We first introduce in Section II the Axion Maxwell equations for a semi-infinite magnetoelectric. In order to better understand the nature of the actual topological defect being induced by the point charge on the TI Section III we derive an exact solution for a vortex of finite length L in the London limit and obtain the Pearl vortex solution in the limit L → 0 and large planar distances compared to the London penetration depth. Such a Pearl vortex regime will arise as a topological defect induced by a point charge sitting precisely on the TI surface. This result will be obtained as a special case of the more general solution in Section IV. In Section III C we compare the solution of a very thin solenoid of finite length L to the obtained solution for a vortex of finite length and notice crucial differences. At large distances the field profiles arising from the end points of the solenoid behave precisely as magnetic monopoles. In the vortex case the monopole-like field profiles do not quite resemble actual monopoles due to Meissner screening inside the superconductor. Such a screening is absent in the case of a solenoid. The general solution of the electric and magnetic field and the associate angular momentum for a charge placed at a distance d away from such a semi-infinite magnetoelectric medium is detailed in Section IV. In Section V B we derive the zero energy mode solutions of the Dirac equation on the TI surface in the presence of the electromagnetic field induced by the external point charge. A general argument for absence of magnetic monopoles is provided in VI and we end with a brief summary and conclusions.

II. AXION MAXWELL EQUATIONS FOR A SEMI-INFINITE MAGNETOELECTRIC
Given the textbook nature of the problem on one the hand and the importance of its exact solution on the other, we present the steps to obtain a direct solution for the electric and magnetic fields in a semi-infinite threedimensional Maxwellian magnetoelectric in some detail. The effective Lagrangian density is given by [1], where L a is given above and Gaussian units are being used. Similarly to Ref. [12], we assume that the (toplogical) magnetoelectric medium occupies the region z ≤ 0, with the surface at z = 0 separating it from a trivial insulator, which we assume to be the vacuum, see Fig. 1. Thus, we have a dielectric constant = 1 and θ = 0 for z > 0. We further assume for simplicity that the magnetic properties are such that µ = 1 for all z. The easiest way to obtain the field equations is to write the standard Maxwell equations in the presence of matter and recall that, We obtain in this way the general field equations in the form, while the source-free Maxwell equations remain unchanged, since the latter are actually constraints following from the Bianchi identity for the electromagnetic field strength.
Die to the second Eq. (2) above, H = B + (αθ/π)E, we see that one way of solving the problem shown in Fig. 1 is to use a magnetic image charge [12], since the problem becomes analogous to the one of a semi-infinite system with a magnetic permeability µ = 1 [17]. Indeed, here the role of magnetization is played by −αθE/(4π 2 ), so we can use the constraint ∇ · B = 0 to define a magnetic charge density, ρ M = ∇ · [αθE/(4π 2 )], leading to the equations ∇ · H = 4πρ M and ∇ × H = 0, which is formally identical to a problem in electrostatics [17,18]. We will not follow this approach here and proceed to solve the equations using a vector potential, as this will help us to clarify the similarities and differences from these image magnetic charges and actual magnetic monopoles.
After setting E = −∇φ and B = ∇×A, we obtain the differential equations for the scalar and vector potentials, where we have assumed the Coulomb gauge ∇ · A = 0. Since ∇θ = −θδ(z)ẑ for the system under consideration, we have to actually solve the equations, subjected to boundary conditions reflecting the discontinuities in the normal derivatives of the potentials at z = 0, implied also by the change of θ at the interface. Translational invariance in the xy-plane implies, whereÂ(p, z) is the Fourier transform of the vector potential in the plane. The above equation is to be solved with the boundary conditions, dÂ dz (12) where η → 0+. One finds after a straightforward calculation that which yields At this point it is important to observe that the magnetic field is divergence-free everywhere in space as the above expression obviously satisfies ∇ · B = 0 everywhere, irrespective of the form of the electric field.

III. VORTEX OF FINITE LENGTH, SOLENOIDS AND ARTIFICIAL MONOPOLES
In order to put better in perspective what kind of solution is actually obtained and understand how objects that from large distance may look like a magnetic monopole, we first consider a vortex line of finite length in the London limit. In the limit where the length of the vortex line approaches zero, one obtains the field of a so called Pearl vortex at large planar distances, r λ L , where λ L is the London penetration depth. We will see later that for a TI a Pearl vortex-like solution [15] is obtained for the case where the point charge sits precisely at the TI surface. The major difference to the actual Pearl vortex in a thin superconducting slab is that the solution for the TI will be exact and not only valid at large distances in the plane.

A. Brief review of London theory
In the London limit the static superconducting current is given by the well known formula [19], where ρ s is the superfluid density and v s is the superfluid velocity, which in the case of a superconductor is given by, where ϕ is the phase of the superconducting order field. For later use, we are labeling the planar coordinates r = (x, y) separately from z. Thus, the London equation is simply the Maxwell equation, supplemented by the constraint, ∇ · B = 0. If the space is simply connected, we have, ∇ × ∇ϕ = 0 everywhere, and the London equation simplifies to, where, The latter equation yields the London penetration depth λ L = m −1 L . Vortices make space multiply connected and the London equation becomes more interesting. Generally a superconductor features both open and closed (loops) vortex lines [20]. If vortex lines are accounted for, the phase gradient has the form [21], where R = (r, z). In the above equation ϕ L denotes the longitudinal part of the phase satisfying ∇ × ∇ϕ L = 0, while the second term corresponds to the contribution from vortex lines L i , 1 ≤ i ≤ N v , with the integral being along the i-th vortex line determined by the vector γ i . Due to the second term, we have ∇ × ∇ϕ = 0. For the simple case of a single infinite straight vortex line, we have γ(z ) = z ẑ, z ∈ (−∞, ∞), such that we obtain, Thus, where n ∈ Z is the vorticity winding number. From Stokes theorem we can therefore write, and the London equation becomes [19], where Φ 0 = hc/(2e) is the elementary flux quantum for a superconductor, corresponding to half of the usual flux quantum arising in the Aharonov-Bohm effect. Therefore, for an infinite system the ANO solution [22,23] for a single infinite vortex line is given in the London limit by [19], where K 0 (x) is a modified Bessel function of second kind.

B. London theory in a superconducting slab
For a superconducting slab of thickness L general vortex line solutions in the London regime have been obtained by Brandt [24] and Carneiro and Brandt [25]. For a straight vortex line parallel to the z-axis we have that vector potential now depends on z and has the form, We assume that the slab occupies the region Thus, the London equation has to be solved in cylindrical coordinates using boundary conditions for A(r, z) both at z = 0 and at z = −L. These are continuity of A(r, z) of its derivative with respect to z at the surfaces z = 0 and z = −L.
The region outside the slab is assumed to be vacuum. The easiest way to obtain the solution is to recall the vector potential leading to the magnetic field (67) and generalize it to, where J 1 (x) is a Bessel function and a(p, z) is determined by the boundary conditions. The infinite system has a solution corresponding to a(p, z) = 1. In this case the integral can be performed exactly to obtain, which yields the vector potential whose curl produces the magnetic field (67).
For a system of thickness L we obtain on the other hand, A regime particularly interesting for us is the large distance one, where r λ L . In this regime we obtain, and thus, corresponding to the magnetic field of a magnetic monopole of charge g = nΦ 0 /(2π) located at R 0 = 0 for z > 0 and to the magnetic field of a magnetic monopole of charge −g located at R 1 = (0, 0, −L) for z < −L. If in addition we consider the thin film limit L → 0 we obtain a point vortex field profile, a so called Pearl vortex [15], Note that despite exhibiting for both z > 0 and z < −L fields of magnetic monopoles, the interval restrictions guarantee that ∇ · B = 0, as it should.
Another limit case of interest is the short distance one corresponding to r λ L . In this case, which leads to the typical thin solenoid expression for the magnetic field in the slab region R, C. Comparison with a thin solenoid of length L The asymptotic behavior of the vortex solution for a slab geometry exhibited some similarities with a solenoid. There are some crucial differences, however, which we analyze below.
It is a well known fact in classical electrodynamics that a thin solenoid can be thought as a line of point magnetic dipoles, an approximation valid at large distances [17]. The vector potential is assumed to be a sum over infinitesimal elements of vector potential associated to an element of magnetic moment, dm = m 0 dz ẑ, where m 0 is the magnetic moment per unit length of the solenoid, assumed to lie along the z-axis. Thus, The vector potential for a thin solenoid of length L is therefore given by integrating the above equation in z ∈ (−L, 0), This yields the magnetic field, where θ H (x) is the Heaviside unit step function. Thus, a thin solenoid of length L at large distances appears as two magnetic monopoles of charges ±m 0 connected by a Dirac string of length L. Note that the presence of the string guarantees that ∇ · B = 0 [26]. The magnetic field (37) is clearly very different from the magnetic field of a vortex line at large distances, Eq. (31). In the latter equation the magnetic fields of the monopoles are completely screened inside the slab. This is not the case in Eq. (37). Also the limit L → 0 does not yield the magnetic field of a Pearl vortex, yielding instead, Interestingly, the result for a straight thin solenoid can be generalized to any curved shape. This can be done by exploring the similarity of Eq. (35) with Eq. (20) for the phase gradient in a superfluid. Accordingly, we consider a line of infinitesimal magnetic moments along a curved line C defined by the vector function, γ(s), s ∈ [0, 1], i.e., dm = m 0 dγ. Thus, instead of Eq. (35) we have, so that, The magnetic field is therefore derived as follows, which immediately leads to, where γ 2 = γ(1) and γ 1 = γ(0). Equation (42) obviously includes Eq. (37) as a special case. Indeed, in this case we have simply γ(z ) = z ẑ, with γ 2 = 0 and γ 1 = −Lẑ, while for the delta function term we have,

A. Calculation of the electric field
We now consider for the charge density a point charge q at z = d > 0 as indicated in Fig.1. Since the Poisson equation is translation invariant in the xy-plane, one can perform a two-dimensional Fourier transform to obtain the differential equations for the Fourier-transformed potential,φ(p, z), The equations for the electric potential Eqs. (44,45) have to obey the four boundary conditions where κ = (1/2)(αθ/π) 2 . The boundary condition 2 above follows directly by inserting Eq. (14) into the Poisson equations, Eqs. (6) and (7) and performing a Fourier transform in the plane. These boundary conditions are used to determine the unknown coefficients by matching the solutions in three regions, After determining A, B, C, and D, we obtain, where z 0 ∈ R, we easily obtain the electric potential, yielding in turn the electric fields for z < 0 and z > 0, Unremarkably, the above expressions reduce to the standard textbook ones when θ = 0. Note that only E(r, z = 0) is needed to determine the magnetic field via Eq. (14).
In view of the axion term, the electric field is discontinuous at z = 0, as evidenced by the boundary conditions above. Thus, we have, Therefore,

B. Calculation of the magnetic field
The most straightforward way to calculate the magnetic field is by inserting the expression for the electric field directly in Eq. (27) and performing the resulting integral. The calculation is considerably easier using a Fourier transform on the TI surface, in which case we obtain from Eq. (27), where we have assumed q = N e, N ∈ Z, and, with Φ 0 = hc/e being the elementary flux quantum. Thus, we obtain the Fourier representation for the vector potential, where, Performing the angular integration in the above equation yields, where J 0 (x) is a Bessel function. Thus, Thus, by performing the inverse Fourier transform, we obtain, whose curl yields, where H is the Heaviside step function. It immediately follows that, since H(−d) = 0. Although the total magnetic flux through any closed surface containing a sphere of radius d centered at the origin vanishes, we note that the flux through the TI surface yields precisely Φ for all z.
For d → 0 we obtain from Eq. (66) that We see that if |z| were replaced by z in Eq. (69), it would precisely yield the vector potential of a straight vortex where for a superconductor Φ = nΦ 0 /2. An important difference between the magnetic field above and the actual Pearl vortex arising in superconductors is that the former holds for all r, no matter small, while the actual Pearl vortex field profile of Eq. (32) follows from the large distance limit r λ L for a flux line of vanishing length. In other words, Eq. (70) describes an exact Pearl vortex.
By removing the sgn(z) factor in Eq. (70) we obtain precisely the magnetic field of a Dirac magnetic monopole: in other words the magnetic field (70) behaves as a monopole for z > 0 and as an anti-monopole for z < 0, yielding in this way ∇ · B = 0, see Fig.2.
The stream density plot associated to the magnetic field components above is shown in Fig. 3 for the reduced coordinates z/d and r/d. We note the presence of an extended solitonic object near z = 0, indicating that the point vortex becomes for d = 0 a kind of pancake vortex.
The magnetic field (67) can ibe interpreted as corresponding to image magnetic charges of strength g ± = ±N Φ/(2π) located at z ± = ∓d, respectively. Thus., the magnetic charge g + at z = −d mirrors the magnetic field at z > 0, while the magnetic charge g − at z = d mirrors the magnetic field at z < 0. However, we have seen in Eq. (68) that ∇ · B vanishes everywhere.

V. ANGULAR MOMENTUM
A. Electromagnetic field contribution to the angular momentum The electromagnetic momentum density Π = (E × B)/(4πc) is nonzero for all z > 0, vanishing otherwise. Thus, the angular momentum carried by the electromagnetic field is obtained by integrating over all space the angular momentum density, (r + zẑ) × Π. The components L x and L y vanish due rotational invariance in the plane, so we obtain, Performing the integral above yields, As with angular momentum of dyons [16], the above angular momentum is independent of the distance of the point charge to the TI surface. However, in contrast to the dyon, the above angular momentum is neither an integer nor a half-integer multiple of . In order for this to happen, Φ = nΦ 0 , n ∈ Z, would have to be satisfied. But this would then in general imply a negative dielectric constant for all n as according to Eq. (61). Interestingly however, L z features a quantization ∼ N 2 characteristic of angular momentum of Chern-Simons vortices [27]. Despite this similarity, it is important to emphasize that the integer squared arising in the latter corresponds to the vortex quantum number rather than the number of particles N as used here in connection with the point charge q = N e. This fact has experimental consequences and will come back to this point below.
It might be tempting to associate the result of Eq. (72) to a dyon exhibiting anyon behavior [12]. As Φ ∝ (1 + + κ) −1 with the dielectric constant of the magnetoelectric medium and κ ∝ α 2 this would imply a dependence of L z on the dielectric and fine structure constant. Such non-universal behavior in the process of exchanging two particles is rather not expected in quantum statistics.
Let us compare the result of Eq. (72) with the one obtained for a vortex at the interface of a SC-TI heterostructure. For the latter it has been demonstrated that [9], where n is the vortex quantum number. The integer number square appearing in both Eqs. (72) and (73) have a completely different origin. In Eq. (72) N refers to the number of electron charge units contained in the charge q at a distance d from the TI surface, while n in Eq. (73) is a winding number associated to flux quantization in superconductors. It turns out that n is rarely larger than the unity. On the other hand, it is easily possible to have considerably larger values of N . Thus, in spite of the small factor α 2 in Eq. (73), L z is not necessarily small and it can actually be of the same order or even larger than the induced angular momentum in SC-TI systems. Another important difference between Eq. (73) and Eq. (72) is that the former refers to the total angular momentum rather than just the one due to the electromagnetic field. Equation (73) actually contains also the mechanical contribution, i.e., one corresponding to the orbital motion of Cooper pairs around the vortex. In the next subsection we will address the contribution of the surface fermions to the total angular momentum.
Interestingly, conservation of angular momentum would imply that the induced magnetoelectric angular momentum be counterbalanced by a mechanical rotation of the TI. This can be in principle be detected experimentally.

B. Angular momentum contribution from Surface Dirac fermions
The Dirac equation on the TI surface is given by, with the Hamiltonian, where v F is the Fermi velocity and in terms of Pauli matrices α = (−σ y , σ x ). The Pauli matrix σ z multiplies ∆ = −µ B L z (µ B is the Bohr magneton), which is a Zeeman term induced by the angular momentum of the electromagnetic field. Were we to have φ = 0 and ∆ = 0, the state σ z ψ E would have been an eigenstate of the Dirac operator with energy −E, and therefore σ z ψ E = ψ −E [28]. This property fails when φ = 0, which is precisely the case we are dealing with here. This implies that any underlying zero modes in this problem are not selfconjugate.
There are two distinct situations to be considered here, depending on whether d = 0 or not. From Eq. (66), we have, so we see that for d = 0 the vector potential vanishes for r → 0 and behaves as A(r, 0) ∼φ/r for r d. For d = 0, on the other hand, the potential is singular at r = 0, but the Dirac equation can still be solved by employing a singular gauge transformation.
Let us consider the case where the charge lies exactly on the TI surface, so that we can set d = 0. In this case, performing the unitary transformation U = e iN (Φ/Φ0)ϕ casts the Hamiltonian in the form, and induces a multivalued phase factor in the spinor field, ψ E = e iN (Φ/Φ0)ϕ Ψ E . Hence, the Dirac equation becomes, , where n is an integer. Note that Ψ E satisfies anti-periodic boundary conditions. It follows that ψ E is an eigenstate of the angular momentum operator, and thus, The total angular momentum of the system is given by adding the above result to the electromagnetic field contribution of Eq. (72), The Ansatz (79) leads to the pair of equations, where we have defined Ω = ω/v F and β = N bΦ/Φ 0 = 2N (c/v F )α/(1+ +κ). Because of the way the equations are coupled, it is more convenient to find solutions by introducing the new variables F = f + g and G = f − g, in which case the equations are recast in the form, Normalized zero mode solutions (Ω = 0) are obtained in the form F (r) = F 0 e −mr r s and G(r) = G 0 e −mr r s , assuming m > 0. These normalized solutions are obtained for s = n 2 − β 2 − 1/2. The boundary conditions at r = 0 require that The constants f 0 = (F 0 + G 0 )/2 and g 0 = (F 0 − G 0 )/2 are easily determined from the normalization condition to be given by, Finally, let us comment on an aspect of zero modes that is usually true in most topological systems but that does not hold in the present case, namely, the existence of an index theorem [29] or, more specifically, the Atiyah-Patodi-Singer (APS) index theorem [30,31]. This theorem states that the number of positive energy modes minus the negative energy ones is an integer topological invariant, the so called η-invariant. In order for this theorem be applicable it is necessary to be able to map oneto-one the positive energy modes to the negative ones, with the zero modes being self-conjugate. When φ = 0 this mapping is provided by the Pauli matrix σ z , similarly to the situations encountered in Refs. [28] and [32] (see also Ref. [29]), where the Pauli matrix σ z acts as conjugation matrix. As already mentioned, in our case σ z does not map ψ E into ψ −E , so we are unable to apply the APS index theorem in this case.

VI. GENERAL ARGUMENT FOR ABSENCE OF INDUCED MAGNETIC MONOPOLES
The constraint ∇ · B = 0 not being violated in the general axion magnetoelectric screening problem can also be argued on the basis of more general considerations. Due to the ∇ · B = 0 constraint a string singularity has to be attached to a Dirac monopole. Monopoles without strings are only possible if topologically nontrivial gauge transformations are allowed [33], in which case two nonsingular vector potentials can be used, A ± (r, z) = ±gr −2 (ẑ × r)(1 ∓ z/ √ r 2 + z 2 ), defined in the regions of a sphere around a point monopole g excluding the south and north poles, respectively. These gauge potentials differ by a singular gauge transformation, since A + − A − = 2g∇ϕ, ϕ ∈ [0, 2π]. Indeed, ∇ × ∇ϕ = 2πδ 2 (r)ẑ. As a consequence, ∇ · B = 4πgδ 2 (r)δ(z). The corresponding topologically nontrivial gauge transformation is therefore G = exp[2iegϕ/( c)], which leads to the Dirac condition, eg/( c) = n/2, n ∈ Z. Such a scenario is not realizable within the axion electrodynamics discussed here, where gauge transformations are topologically trivial. In this case magnetic monopoles would require non-vanishing currents at large distances, contradicting one of the basic tenets of electromagnetism.

VII. CONCLUSION
We have obtained the induced magnetic field due to a charged particle above the surface of a topological insulator or any magnetoelectric material in general. The exact magnetic field has been obtained directly without using image magnetic charges. The solution allows nevertheless for a clear identification of the reflection and transmission image magnetic charges. However, the latter cannot be interpreted as induced magnetic monopoles, since the vector potential does not allow for topologically nontrivial gauge transformations and there is no flux tube connecting the magnetic image charges. In the limit case where the charge lies exactly at the surface, the field of a point vortex, also known as Pearl vortex [15] is obtained. Such a point vortex resembles a monopole, but it is quite different from it, as it does not allow for topologically nontrivial gauge transformations necessary to make ∇ · B = 0.
We have found that the vortex-like solution features a nontrivial angular momentum for the electromagnetic field, in a situation reminiscent to the vortex solution for a TI proximate to a type II superconductor obtained in Ref. [9]. As in the latter reference, we calculated the angular momentum exactly. However, there is an essential difference between them, namely, the result we have obtained in Eq. (72) has a dependence on the dielectric constant, while such a dependence cancels out in the calculation leading to Eq. (73), corresponding to the results obtained in Ref. [9]. The lack of dependence on the dielectric constant in Eq. (73) is due to the fact that this angular momentum corresponds to both mechanical and electromagnetic contribution. However, accounting for the mechanical contribution from the surface Dirac fermions in the presence of an external point charge does remove the -dependence in the total angular momentum. For this reason, one cannot claim that the system exhibits fractional statistics, since quantum statistics should be a universal property of the system and not be dependent on specific material properties like the dielectric constant. Nevertheless the result of an induced total angular momentum is a very interesting one, since it indicates that a point charge in proximity to a TI would induce a compensating mechanical torque that can be detected experimentally.
Since the Dirac equation is known to favor zero mode solutions in the presence of topological defects like vortices or monopoles [32], we have also investigated this possibility here and obtained that zero modes indeed exist. However, the fermions have to be gapped and thus break TR, otherwise no normalized zero modes can be found.