# Shadows and lensing of black holes immersed in strong magnetic fields

Haroldo C. D. Lima Junior<EMAIL_ADDRESS>Faculdade de Física,
Universidade Federal do Pará, 66075-110, Belém, PA, Brazil  Pedro V. P. Cunha
<EMAIL_ADDRESS>Departamento de Matemática da Universidade de Aveiro and Centre
for Research and Development in Mathematics and Applications (CIDMA), Campus
de Santiago, 3810-183 Aveiro, Portugal  Carlos A. R. Herdeiro<EMAIL_ADDRESS>Departamento de Matemática da Universidade de Aveiro and Centre for Research
and Development in Mathematics and Applications (CIDMA), Campus de Santiago,
3810-183 Aveiro, Portugal  Luís C. B. Crispino<EMAIL_ADDRESS>Faculdade de
Física, Universidade Federal do Pará, 66075-110, Belém, PA, Brazil

###### Abstract

We investigate the null geodesic flow and in particular the light rings (LRs),
fundamental photon orbits (FPOs) and shadows of a black hole (BH) immersed in
a strong, uniform magnetic field, described by the Schwarzschilld-Melvin
electrovacuum solution. The empty Melvin magnetic Universe contains a tube of
planar LRs. Including a BH, for weak magnetic fields, the shadow becomes
oblate, whereas the intrinsic horizon geometry becomes prolate. For strong
magnetic fields (overcritical solutions), there are no LRs outside the BH
horizon, a result explained using topological arguments. This feature,
together with the light confining structure of the Melvin universe yields
panoramic shadows, seen (almost) all around the equator of the observer’s sky.
Despite the lack of LRs, there are FPOs, including polar planar ones, which
define the shadow edge. We also observe and discuss chaotic lensing, including
in the empty Melvin universe, and multiple disconnected shadows.

##  I Introduction

The bending of light rays due to a gravitational field is one of the most
important predictions of general relativity (GR). For a body with the Sun’s
mass, the deflection, first predicted by Einstein, is less than two
arcseconds. After some frustrated attempts [1], the experimental measurement
of the bending of a light ray passing close to the Sun was first accomplished
during the 1919 total solar eclipse. The deviation angle was measured to be
approximately 1.75 arcseconds, and it was in agreement with the values
predicted by GR [2] (cf. Ref. [3] for a historical account). Higher deflection
angles (around ten arcseconds) are obtained for the light emitted by quasars
(active galactic nuclei with high luminosity, discovered during the 1960s) and
deflected by galaxies along the line of sight.

The deflection angles for the Sun and quasars are relatively small. For black
holes (BHs), these angles can be arbitrarily large, and light rays can even
describe planar circular closed orbits dubbed in the literature as light rings
(LRs). It has been proven that stationary, axisymmetric and asymptotically
flat BHs must have one “standard" LR outside the event horizon for each
rotation sense [4].111See also [5; 6; 7]. Such LRs are directly related to the
shadows of BHs [8]. The shadow of a BH [9] is the dark region formed when a BH
is illuminated by some light source, for instance by a celestial sphere, as
seen by an observer inside this celestial sphere. The first proper analysis of
the shadow cast by an isolated spinning BH (described by the Kerr metric [10])
was performed by Bardeen [11]. Over the last few years, in particular
motivated by the Event Horizon Telescope (EHT) collaboration imaging of the
shadow of M87* [12], the computation of BH shadows has become an active area
of research, as the means to test GR and modified gravity, the BH paradigm,
and to also gain insight into strong gravity features and into new physics,
see, $e.g.$, for some recent work, Refs. [13; 14; 15; 16; 17; 18; 19; 20; 21;
22; 44; 45; 46; 23; 24; 25; 26; 27; 28; 29; 30; 31; 32; 33; 34; 35; 36; 37;
38; 39; 40; 41; 42; 43].

In astrophysically realistic environments, BHs are not isolated. Instead, they
are surrounded by matter and fields that can interact with the central BH,
giving rise to interesting phenomena. For instance, the magnetic field
generated due to matter accretion by the BH allows energy extraction through
the so-called Blandford-Znajek mechanism [47]. Moreover, BHs may be subjected
to an external magnetic field coming from a neutron star companion close to
the BH. This has been observed to be the case of Sagittarius A* (Sgr A*), in
the center of the Milky Way, which is surrounded by a magnetar known as SGR
J1745-29 [48; 49; 50]. The magnetic field of SGR J1745-29 is $B=1.6\,\times
10^{14}$ gauss, being among the strongest magnetic fields observed so far
[51]. A recently released EHT analysis, moreover, corroborates the existence
of strong magnetic fields near the event horizon of M87* [52].

These environments, where an external strong magnetic field surrounds a BH,
can be modeled analytically with the Schwarzschild-Melvin BH (SMBH) exact
solution of Einstein-Maxwell theory, obtained by Ernst [53]. The SMBH solution
represents a Schwarzschild BH spacetime immersed in a Melvin universe [54].
The latter is a self gravitating “constant" magnetic field and it is not
asymptotically flat, since the magnetic field extends all the way to spatial
infinity. Despite this unphysical feature, in what concerns astrophysical BHs,
the SMBH geometry can, nonetheless, be relevant in describing the near horizon
geometry of a BH subjected to certain strong magnetic fields.

The SMBH spacetime presents interesting properties, which have been
extensively studied. For instance, the motion of charged particles on the
equatorial plane of the SMBH spacetime was studied in Refs. [55; 56; 57; 58],
while the motion of timelike and null geodesics on the equatorial plane was
analyzed in Refs. [59; 60; 61], and it was noticed that the SMBH spacetime can
have no LR at all outside the event horizon. One of the goals of this paper is
to explore this intriguing result from the viewpoint of the topological
arguments provided in [4], that guarantee that asymptotically flat BHs must
have LRs. Another goal is to understand how this result impacts on the shadow
of these BHs.

Rotating Kerr BHs in the presence of a magnetic field were studied by Wald
[62], where the rotation axis is aligned to the magnetic field. Interestingly,
this rotating Kerr BH surrounded by an external magnetic field acquires a net
electric charge $Q=2BJ$, where $B$ is the strength of the magnetic field and
$J$ is the angular momentum of the BH. Moreover the absorption and scattering
of planar scalar waves in the SMBH in the limit of weak magnetic field [63],
the quasinormal modes of the static and rotating BHs in external magnetic
fields [64; 65; 66] as well as stationary scalar clouds [67] were also
studied.

As the mass parameter of the SMBH goes to zero, the spacetime reduces to the
Melvin magnetic Universe [54], which was originally proposed as an example of
pure electromagnetic geon, i.e. an electromagnetic field held together by its
own gravitational field [68]. The properties of the motion of charged and
uncharged massive particles, as well as photons, in the Melvin spacetime were
also investigated in [69; 70; 71], which was generalized to Einstein-nonlinear
electrodynamics in [72].

In this paper we shall analyze the LRs and shadows of static chargeless BHs
surrounded by a strong external magnetic field, modeled by the SMBH spacetime.
As a particular (vanishing mass parameter) case, we study the LRs and the
gravitational lensing of the Melvin magnetic Universe. The remainder of this
paper is organized as follows. In Sec. II we review the properties of the SMBH
spacetime. In Sec. III we define a 2D effective potential, study the motion of
null geodesics and the LRs in the SMBH spacetime, as well as in the Melvin
solution. In Sec. IV we study the shadows cast by the SMBH, the gravitational
lensing of the Melvin spacetime and present a selection of our backwards ray-
tracing results. Our final remarks are presented in Sec. V.

##  II The spacetime

The geometry of the SMBH spacetime, which is a static and axially symmetric BH
solution of the Einstein-Maxwell equations, is described by the following line
element [53]:

$\displaystyle ds^{2}=-$
$\displaystyle\Lambda^{2}(r,\theta)\left(1-\frac{2M}{r}\right)\,dt^{2}+\frac{\Lambda^{2}(r,\theta)}{\left(1-\frac{2M}{r}\right)}\,dr^{2}$
$\displaystyle+$ $\displaystyle
r^{2}\Lambda^{2}(r,\theta)\,d\theta^{2}+\frac{r^{2}\,\sin^{2}\theta}{\Lambda^{2}(r,\theta)}\,d\phi^{2},$
(1)

where

$\displaystyle\Lambda(r,\theta)=1+\frac{B^{2}r^{2}\sin^{2}\theta}{4},$ (2)

$M$ is the mass of the BH and $B$ determines the strength of the external
magnetic field.222In this paper we assume that $B\geqslant 0$. This choice
does not imply in any loss of generality. The electromagnetic potential 1-form
is given by [58]

$A=\frac{Br^{2}\sin^{2}\theta}{2\Lambda}d\phi.$ (3)

Due to the external “constant" magnetic field, the SMBH geometry with $B\neq
0$ does not approach the Minkowski spacetime in the asymptotic region. For
$B=0$, the SMBH spacetime reduces to the (vacuum) Schwarzschild solution,
while for $M=0$ it reduces to the Melvin magnetic Universe [54].

The line element (II) exhibits singularities at $r=0$ and $r=2M$. The former
corresponds to an irremovable singularity, since the Kretschmann scalar

$K=\frac{\mathcal{P}(r,\theta)}{r^{6}\Lambda^{8}},$ (4)

diverges at $r=0$, regardless of the value of the angular coordinate $\theta$.
In Eq. (4), $\mathcal{P}(r,\theta)$ is a long and nonenlightening function of
the coordinates $(r,\theta)$, obeying

$\lim_{r\rightarrow 0}\mathcal{P}(r,\theta)=48M^{2}.$ (5)

The singularity at $r=2M$ is merely a coordinate one. This can be confirmed by
writing the SMBH line element (II) in terms of the Eddington-Finkelstein-like
coordinate $u$, defined through

$du=dt+\frac{dr}{1-\frac{2M}{r}},$ (6)

obtaining that

$\displaystyle ds^{2}=$
$\displaystyle-\Lambda^{2}\left(1-\frac{2M}{r}\right)du^{2}+2\Lambda^{2}dr\,du+r^{2}\Lambda^{2}d\theta^{2}$
$\displaystyle+\frac{r^{2}\,\sin^{2}\theta}{\Lambda^{2}}\,d\phi^{2},$ (7)

which is regular at $r=2M$. In the Schwarzschild geometry, $r=2M$ corresponds
to the location of the BH event horizon. The event horizon for asymptotically
flat geometries is a well-defined teleological concept: the boundary of the
causal past of future null infinity. Identifying it requires knowledge of the
global structure of the spacetime [73]. However, the formulation for
nonasymptotically flat spacetimes may introduce conceptual and technical
difficulties. The existence of a BH region in the SMBH spacetime may,
nonetheless, be justified by the identification of an apparent horizon, which
requires only local knowledge [73]. To show $r=2M$ delimits a BH region, in
this sense, let us denote the tangent vectors of a congruence of outgoing and
ingoing null geodesics by $l^{\mu}$ and $n^{\mu}$, respectively. These
geodesic congruences are not necessarily affinely parametrized, i.e.,

$\displaystyle l^{\mu}\nabla_{\mu}l^{\nu}=k_{l}l^{\nu},$ (8) $\displaystyle
n^{\mu}\nabla_{\mu}n^{\nu}=k_{n}n^{\nu},$ (9)

where $k_{l},k_{n}$ are some functions. The expansion $\Theta$ of these null
congruences are given by

$\displaystyle\Theta_{(l)}=\nabla_{\mu}l^{\mu}-k_{l},$ (10)
$\displaystyle\Theta_{(n)}=\nabla_{\mu}n^{\mu}-k_{n}.$ (11)

The apparent horizon is defined as the outer boundary of a trapped region,
where the following condition holds [73; 74; 75]:

$\displaystyle\Theta_{(l)}=0.$ (12)

For radially outgoing and ingoing null geodesic congruences along the
equatorial plane, we have

$\displaystyle
l^{\mu}=\left(l^{u},l^{r},l^{\theta},l^{\phi}\right)=\left(2,1-\frac{2M}{r},0,0\right),$
(13) $\displaystyle
n^{\mu}=\left(n^{u},n^{r},n^{\theta},n^{\phi}\right)=\left(0,-2,0,0\right),$
(14)

respectively. Using Eqs. (13) and (14) together with Eqs. (8) and (9), we find
that

$\displaystyle
k_{l}=\frac{2M}{r^{2}}+\frac{B^{2}\sin\theta^{2}}{\Lambda}\left(r-2M\right),$
(15) $\displaystyle k_{n}=-\frac{2Br\sin\theta^{2}}{\Lambda}.$ (16)

Moreover, from Eqs. (10), (13), and (15), we have

$\displaystyle\Theta_{(l)}=\frac{2\left(r-2M\right)}{r^{2}},$ (17)

which vanishes at $r=2M$, as required for an apparent horizon. For the sake of
completeness, we can also show that the expansion of null ingoing geodesics is
negative

$\Theta_{(n)}=-\frac{4}{r}<0.$ (18)

Therefore, we conclude that $r_{h}=2M$ is the radial coordinate of an apparent
horizon in the SMBH spacetime, since at this point the expansion (17) is zero
and for $r<2M$ it is negative.

The surface area of the apparent horizon is

$A_{h}=\left.\int_{0}^{2\pi}\int_{0}^{\pi}\sqrt{g_{\theta\theta}g_{\phi\phi}}d\theta
d\phi\right|_{r=r_{h}}=4\pi r_{h}^{2},$ (19)

which is equal to that of the vacuum Schwarzschild geometry. Although the
radial coordinate and the surface area of the apparent horizon are equal in
the SMBH and vacuum Schwarzschild case, the geometry of the apparent horizon
is different, since the SMBH spacetime is not spherically symmetric. To
illustrate this point, in Fig. 1 we show the isometric embedding in Euclidean
3-space of the SMBH apparent horizon geometry [76; 77; 78]. The apparent
horizon keeps a $\mathbb{Z}_{2}$ north-south symmetry with a well-defined
equator, but becomes increasingly prolate and flattened near the equator, as
$BM$ is increased. In the next sections, we shall explore the properties of
null geodesics on and outside the equatorial plane.

Figure 1: Embedding of the apparent horizon geometry in Euclidean 3-space for
the SMBH spacetime with several values of $B$. For $B=0$, it reduces to the
vacuum Schwarzschild case. For $B>0$, the apparent horizon geometry becomes
increasingly prolate.

##  III Null geodesics

###  III.1 LRs in the SMBH (and Melvin) spacetimes

The motion of null geodesics can be described by Hamilton’s equations:

$\displaystyle\dot{x}^{\mu}=\frac{\partial\mathcal{H}}{\partial p_{\mu}},$
(20) $\displaystyle\dot{p}_{\mu}=-\frac{\partial\mathcal{H}}{\partial
x^{\mu}},$ (21)

where the Hamiltonian is given by

$\mathcal{H}=\frac{1}{2}g^{\mu\nu}p_{\mu}p_{\nu}=0.$ (22)

Since the metric components of the SMBH do not depend on the $t$ and $\phi$
coordinates, we have two conserved quantities, namely

$\displaystyle p_{t}=-E,$ (23) $\displaystyle p_{\phi}=L.$ (24)

The quantities $\\{E,L\\}$ are usually interpreted as the energy and angular
momentum of the light ray. However, due to nontrivial Melvin asymptotics, we
remark that $E$ is not generically equal to the energy measured by an observer
at spatial infinity (nor should such an observer be regarded as preferential).
For instance, the energy measured by a static observer at a radius $r$ and
colatitude $\theta$ is given by

$E_{obs}=-p_{\mu}{u}^{\mu}=\frac{E}{\left(1+\frac{B^{2}r^{2}\sin^{2}\theta}{4}\right)\left(1-\frac{2M}{r}\right)^{\frac{1}{2}}},$
(25)

where ${u}^{\mu}$ is the 4-velocity of the observer. When $B\neq 0$ and the
observer lies outside the symmetry axis, $\theta\neq 0,\pi$,
$E_{obs}\rightarrow 0$ as $r\rightarrow\infty$. If the observer lies along the
axis of symmetry, however, $E_{obs}=E$ at spatial infinity. Hence, one may
describe $E$ as the energy of the photon as measured by an observer lying on
the axis of symmetry, i.e., along the direction of the magnetic field lines,
at $r\rightarrow\infty$.

Since the SMBH spacetime is not (generically) spherically symmetric, the
general motion cannot be reduced to the equatorial plane. Moreover the
Hamilton-Jacobi equation for the SMBH is not known to be separable; hence one
is not able (in general) to write the equations of motion as four first order
differential equations, as, for instance, in the Kerr case [79]. The case
$M=0$ (Melvin solution), however, admits separability due to an extra Killing
vector related to the translation symmetry along the $z$ direction [69].

The Hamiltonian $\mathcal{H}$ can be rewritten as a sum of a kinetic term and
a potential term, namely:

$\displaystyle\mathcal{H}=T(r,\theta)+V(r,\theta,E,L),$ (26)

where

$\displaystyle
T(r,\theta)=g^{rr}(p_{r})^{2}+g^{\theta\theta}(p_{\theta})^{2},$ (27)
$\displaystyle
V(r,\theta,E,L)=-\frac{E^{2}}{\Lambda^{2}\left(1-\frac{2M}{r}\right)}+\frac{\Lambda^{2}L^{2}}{r^{2}\sin^{2}\theta}.$
(28)

Since the potential $V(r,\theta,E,L)$ depends on $E$ and $L$, we define a new
effective potential, $H(r,\theta)$, which is independent of $E$ and $L$, given
by

$H(r,\theta)\equiv\frac{\Lambda^{2}\left(1-\frac{2M}{r}\right)^{\frac{1}{2}}}{r\sin\theta},$
(29)

so that we can write

$V(r,\theta,E,L)=\frac{L^{2}}{\Lambda^{2}\left(1-\frac{2M}{r}\right)}\left(H(r,\theta)+\frac{1}{\eta}\right)\left(H(r,\theta)-\frac{1}{\eta}\right),\\\
$ (30)

where $\eta\equiv L/E$. Since $T\geqslant 0$ and $V\leqslant 0$, the null
geodesics in the SMBH spacetime obey the following inequality:

$\frac{1}{|\eta|}\geqslant H(r,\theta),$ (31)

where the equality holds at turning points. When $B\neq 0$, as
$r\rightarrow\infty$, the effective potential diverges (away from the symmetry
axis), $H(r,\theta)\rightarrow\infty$. Thus, Eq. (31) implies that a generic
light ray (not necessarily restricted to the equatorial plane) can only escape
to infinity if $\eta=0$. This property holds for the whole SMBH family, with
$B\neq 0$, including the $M=0$ Melvin universe.

Let us now turn our attention to LRs, i.e., planar circular photon orbits, by
analyzing the properties of the effective potential (29). LRs are defined by
$p_{r}=p_{\theta}=0$ and $\dot{p}_{r}=\dot{p}_{\theta}=0$.333As we shall see
later, in the SMBH spacetime there are planar null (but geometrically
noncircular) orbits which are polar, rather than equatorial - Fig. 10 -, with
$p_{\theta}\neq 0$ e $\eta=0$. These are not LRs and they do not appear as
critical points of $H(r,\theta)$. The tangent vector to these arbits is not a
combination of the Killing vector fields $\partial_{t}$ and $\partial_{\phi}$,
which is the case for LRs [5]. They correspond to the critical points of $V$,
i.e.:

$\displaystyle V=0,$ (32) $\displaystyle\nabla V=0,$ (33)

which may be rewritten in terms of $H(r,\theta)$ as

$\displaystyle H(r,\theta)=\frac{1}{\eta},$ (34) $\displaystyle\nabla
H(r,\theta)=0.$ (35)

Consider first the Melvin spacetime, by choosing $M=0$ in Eq. (29). From Eq.
(35) it is straightforward to identify one equatorial ($\theta=\pi/2$) LR at

$r=\frac{2}{\sqrt{3}B}.$ (36)

This LR, however, is not unique due to the existence of an additional
spacelike Killing vector in the Melvin universe. This can be best seen by
writing the Melvin solution in ($t$, $\rho$, $z$, $\phi$) coordinates, where

$\displaystyle\rho=r\sin\theta,$ (37) $\displaystyle z=r\cos\theta.$ (38)

The line element (II) becomes

$\displaystyle
ds^{2}=\Lambda^{2}(\rho)\left(\,-dt^{2}+d\rho^{2}+dz^{2}\right)+\frac{\rho^{2}}{\Lambda^{2}(\rho)}d\phi^{2},$
(39) $\displaystyle\Lambda(\rho)=1+\frac{B^{2}\rho^{2}}{4},$ (40)

which is independent of $t$, $z$ and $\phi$. The LR (36), in the
($t,\rho,z,\phi$) coordinates, is determined by

$\displaystyle\rho=\frac{2}{\sqrt{3}B},$ (41) $\displaystyle z=0.$ (42)

Since the Melvin universe has a translation symmetry along the $z$ direction,
there is a LR at each $z=$ constant plane with a cylindrical radius given by
$\eqref{rho_LR}$. Hence, we conclude that the Melvin universe admits a tube of
planar LRs. This means that, instead of just a single LR orbit, there exists
an infinite continuous set of LRs. In Fig. 2 (top panel) we show the effective
potential $H(r,\pi/2)$ for the Melvin spacetime. The tube of planar LRs is
stable against radial perturbations.

Considering now the SMBH spacetime, from the LR conditions, Eqs. (34) and
(35), we find that LRs must satisfy both $\theta=\pi/2$ and

$\displaystyle 3B^{2}r^{3}-5MB^{2}r^{2}-4r+12M=0\ .$ (43)

In Fig. 2 (bottom panel), the effective potential $H(r,\pi/2)$ is plotted. It
may present two local extrema, corresponding to one stable and one unstable
LR, which are determined by Eq. (43). The red dots correspond to unstable LRs,
while the black crosses correspond to stable LRs. In contrast to the Melvin
spacetime, the SMBH does not admit a translation symmetry along the $z$
direction, so that the LRs are located on the equatorial plane only. Due to
the shape of the potential, bound equatorial photon orbits are allowed around
the stable LR. Such orbits neither fall into the BH nor escape to infinity. In
Fig. 3, we show an example of such a bound orbit in the SMBH spacetime with
$BM=0.1$.

Figure 2: Top: Effective potential $H(r,\theta)$ for the Melvin spacetime,
evaluated at the equatorial plane. Bottom: Effective potential for the SMBH
spacetime, also evaluated at the equatorial plane. We have chosen different
values of magnetic field strength. The red (black) dots (crosses) correspond
to the equatorial unstable (stable) LRs radial coordinate. For $B>B_{c}$ there
are no LRs (green curve).  Figure 3: Example of bound orbit on the equatorial
plane with $\eta=5.004$ and $BM=0.1$. The orbit is bounded between $r=3.396M$
and $r=20.963M$. The central circle represents the BH apparent horizon.

Bound orbits are also allowed outside the equatorial plane in the SMBH
spacetime, due to the presence of closed pockets on the effective potential
$H(r,\theta)$, see $e.g.$ [17]. The existence of such closed pockets gives
rise to interesting lensing features as we shall see below.

As the value of $B$ is increased (in terms of the scale $M$), stable and
unstable LRs converge in their radial coordinate — see the sequence of
“crosses” in Fig. 2. For a critical value of the magnetic field strength
($B_{c}$) the radial coordinates of the stable and unstable LRs coincide. The
value of $B_{c}$ is determined by Eq. (43) together with

$\partial^{2}_{r}H(r,\pi/2)=0,$ (44)

which has the solution

$\displaystyle B_{c}M=\frac{2}{5}\sqrt{\frac{169-38\sqrt{19}}{15}}\approx
0.189366.$ (45)

For $B>B_{c}$, Eq. (43) possesses no real solution with $r>2M$, and therefore
there are no LRs outside the apparent horizon. This intriguing fact is further
explored in the next section. In the following, SMBHs with $B>B_{c}$ shall be
called overcritical.

###  III.2 Topological charge of LRs in the SMBH and Melvin spacetimes

It was recently shown that a four-dimensional BH spacetime that is stationary,
axi-symmetric, circular and asymptotically flat must have one standard LR
outside the event horizon for each rotation sense [4]. Subsequently, a
generalization of this result for asymptotically de Sitter/anti-de Sitter
spacetimes was also reported, with a similar conclusion [6]. Since the SMBH
solution is asymptotically Melvin, the existence of LRs is not guaranteed by
the above theorems; thus, the intriguing possibility of a BH without LRs
emerges. As seen in the previous section, this is indeed the case for
overcritical SMBHs.

In order to have a deeper understanding on the absence of LRs in a BH
spacetime, it is insightful to analyze how the different (Melvin) asymptotic
structure impacts on the LR topological charge, which is at the heart of the
aforementioned theorems [5; 4]. To do so, let us introduce the vector field
$\textbf{v}=\left(\text{v}_{r},\text{v}_{\theta}\right)$, where

$\displaystyle\text{v}_{i}\equiv\frac{\partial_{i}H}{\sqrt{g_{ii}}},\ \ \ \ \
\ \ \ \ i=\left(r,\theta\right)\ \ \ \ \ \text{(not summed)}.$ (46)

Thus

$\displaystyle\text{v}_{r}=\frac{\sin\theta}{4r^{3}}\left[B^{2}r^{2}\left(3r-5M\right)-4\left(r-3M\right)\csc^{2}\theta\right],$
(47)
$\displaystyle\text{v}_{\theta}=\frac{\sqrt{1-\frac{2M}{r}}\cos\theta}{16r^{2}\Lambda}\left[B^{2}r^{2}\left(3B^{2}r^{2}\sin^{2}\theta+8\right)-16\csc^{2}\theta\right].$
(48)

We write the components of v in terms of its “norm” $v$ and the angle
$\Omega$:

$\displaystyle\text{v}_{r}=\text{v}\,\cos\Omega,$ (49)
$\displaystyle\text{v}_{\theta}=\text{v}\,\sin\Omega,$ (50)

where

$\displaystyle\text{v}^{2}\equiv\partial^{\mu}H\partial_{\mu}H=\text{v}_{r}^{2}+\text{v}_{\theta}^{2}.$
(51)

Figure 4: Top panel: plot of the normalized vector field v on the $(r,\theta)$
space for the SMBH spacetime with $BM=0.1$. Bottom panel: depiction of a
closed (dotted) curve $C$ that encloses two LRs in the SMBH spacetime. The
curve $C$ is made to approach the axis and horizon under appropriate limits
(see text and [4] for details).

As an illustration of this vector field on the SMBH spacetime, we have
represented in Fig. 4 (top panel) the field v on the $(r,\theta)$ space. When
compared to the Schwarzschild spacetime, the main difference occurs in the
asymptotic region (see Fig. 1 of Ref. [4] for the Schwarzschild case). The
Melvin asymptotics make $\text{v}_{r}$ asymptotically positive, whereas it is
negative for asymptotically flat spacetimes. In contrast, the qualitative
behavior of the vector field next to the horizon is similar for the SMBH and
the vacuum Schwarzschild BH. The asymptotic difference will impact decisively
on the topological charge.

To compute the topological charge, consider a piecewise smooth and positive
oriented curve $C$ in $(r,\theta)$ space; it is possible to introduce the
useful integer quantity [4]:

$w=\frac{1}{2\pi}\oint_{C}d\Omega.$ (52)

If one deforms the curve $C$ without it intersecting any LR, the value of $w$
remains unaltered. Due to this, $w$ can be interpreted as a topological charge
associated to $C$. Given a single nondegenerate LR that is completely inside
the curve $C$, then $w$ is the topological charge of that LR [4]. Following
[4], a LR with $w=-1$ is named standard, while a LR with $w=1$ is named exotic
[5; 4]. The LR stability is actually closely connected to this value: the LR
is stable (unstable) if $w=+1$ ($w=-1$), provided that the null energy
condition is satisfied [5], which is the case in electrovacuum.

The total topological charge of the SMBH spacetime ($M\neq 0$) is computed by
choosing a curve $C$ that encloses all of the spacetime outside the BH, under
some appropriate limit [4]. Such a possible curve $C=L_{1}\cup L_{2}\cup
L_{3}\cup L_{4}$ is depicted in Fig. 4 (bottom panel), where $C$ is
illustrated enclosing two LRs, one stable and the other one unstable. The
paths $\\{L_{4},L_{2}\\}$ have constant $\theta-$values,
$\theta=\\{\varepsilon,\,\pi-\varepsilon\\}$ respectively, whereas
$\\{L_{1},L_{3}\\}$ have constant radial coordinates $\\{r_{0},R\\}$,
satisfying $R>r_{0}>2M$. Thus, (52) can be separated into the sum:

$w=\frac{1}{2\pi}\left(\Delta\Omega_{1}+\Delta\Omega_{2}+\Delta\Omega_{3}+\Delta\Omega_{4}\right),$
(53)

where

$\displaystyle\Delta\Omega_{1}=\left.\int_{L_{1}}\frac{d\Omega}{d\theta}d\theta\right|_{r=R},$
(54)
$\displaystyle\Delta\Omega_{2}=\left.\int_{L_{2}}\frac{d\Omega}{dr}dr\right|_{\theta=\pi-\varepsilon},$
(55)
$\displaystyle\Delta\Omega_{3}=\left.\int_{L_{3}}\frac{d\Omega}{d\theta}d\theta\right|_{r=r_{0}},$
(56)
$\displaystyle\Delta\Omega_{4}=\left.\int_{L_{4}}\frac{d\Omega}{dr}dr\right|_{\theta=\varepsilon}.$
(57)

The LR charge outside the BH is then obtained by taking
$\lim_{R\to\infty}\,\lim_{r_{0}\to 2M}\,\lim_{\varepsilon\to 0}w$. Computing
analytically the integrations (54)-(57) for SMBHs proved to be challenging in
the most general case. However, in the large magnetic field limit, $BM\gg 1$
the calculations simplify considerably. In such limit we obtain:

$\displaystyle\text{v}_{r}\simeq
B^{2}\sin\theta\left(\frac{3r-5M}{4r}\right),$ (58)
$\displaystyle\text{v}_{\theta}\simeq\frac{3}{4}B^{2}\cos\theta\sqrt{1-\frac{2M}{r}}.$
(59)

Taking also $\varepsilon\ll 1$, the integrations (54)-(57) yield:

$\displaystyle\Delta\Omega_{1}\simeq-\pi-2\varepsilon\,\mathcal{F}(R)+\quad\mathcal{O}\left(\varepsilon^{3}\right),$
(60)
$\displaystyle\Delta\Omega_{2}\simeq\varepsilon\Big{[}\mathcal{F}(R)-\mathcal{F}(r_{0})\Big{]}+\quad\mathcal{O}\left(\varepsilon^{3}\right),$
(61)
$\displaystyle\Delta\Omega_{3}\simeq\pi+2\varepsilon\,\mathcal{F}(r_{0})+\quad\mathcal{O}\left(\varepsilon^{3}\right),$
(62)
$\displaystyle\Delta\Omega_{4}\simeq\varepsilon\Big{[}\mathcal{F}(R)-\mathcal{F}(r_{0})\Big{]}+\quad\mathcal{O}\left(\varepsilon^{3}\right),$
(63)

where we have used the auxiliary function

$\mathcal{F}(r)=\frac{(5M-3r)}{3\sqrt{r(r-2M)}}.$

The sum (53) then readily gives us
$w\simeq\mathcal{O}\left(\varepsilon^{3}\right)$. Since the subleading term
can be made much smaller than unity, and accounting for the fact that $w$ must
be an integer, we then arrive at the neat result that $w=0$, in the large
magnetic field limit $(BM\gg 1)$.

Regarding the most general case, it is nevertheless straightforward to check
numerically that the integrations (54)-(57) also yield a vanishing topological
LR charge for the SMBH spacetime if $B\neq 0$:

$\quad w_{\textrm{SMBH}}=0\quad\qquad\textrm{\small$(B\neq 0)$}.$ (64)

In contrast, Schwarzschild has $w_{\textrm{Schw}}=-1$. Since the topological
charge $w$ is constant along the SMBH family with different $B\neq 0$, we may
have one of the following two cases:

  * •

(i) undercritical SMBHs ($B<B_{c}$): A BH spacetime with one standard unstable
LR ($w_{1}=-1$) and one stable exotic LR ($w_{2}=+1$), such that
$w=w_{1}+w_{2}=0$.

  * •

(ii) overcritical SMBHs ($B>B_{c}$): A BH with no LRs at all outside the
horizon ($w=0$).

Case (ii) is the only example (known to the authors) of a BH spacetime without
any LR outside a regular horizon. This is possible due to the nontrivial
asymptotic behavior of the SMBH (for $B\neq 0)$. The absence of LRs in the
SMBH spacetime with $B>B_{c}$ was first discussed in [60].

Figure 5: Plot of the normalized vector field
$\textbf{v}=(\text{v}_{\rho},\text{v}_{z})$ in the $(\rho,z)$ plane, for the
Melvin spacetime with $B=0.1$. The vertical red line represents the location
of the tube of planar LRs.

As a final remark concerning the analysis in this section, the contour
integration approach to compute $w$ for the entire spacetime is not, strictly
speaking, well defined for the Melvin $(M=0)$ case, since any closed curve $C$
setup to include virtually all of the Melvin spacetime must intersect a LR in
the tube of planar LRs. To illustrate this, in Fig. 5 we represent the vector
field v in the $(\rho,z)$ plane for the Melvin spacetime, described by

$\displaystyle\text{v}_{\rho}=\frac{\partial_{\rho}H(\rho,z)}{\sqrt{g_{\rho\rho}}}=\frac{3B^{2}\rho^{2}-4}{4\rho^{2}},$
(65)
$\displaystyle\text{v}_{z}=\frac{\partial_{z}H(\rho,z)}{\sqrt{g_{zz}}}=0.$
(66)

In that plot, the vertical red line represents the tube of planar LRs, which
extends indefinitely along the $z$ direction. This tube creates an obstacle
for the contour setup depicted in Fig. 4 (bottom). Although it might be
possible to still define a topological charge for the Melvin universe, it will
require a different approach which is beyond the scope of this paper.

##  IV Shadows and gravitational lensing

###  IV.1 Setup

Let us now turn to the gravitational lensing in the Ernst spacetime and the
shadows cast by SMBHs. We use the backwards ray-tracing method, which consists
in evolving the light rays from the observer’s position, and backwards in
time, until it is captured by the BH or scattered to a celestial sphere of
finite radius. For this purpose, we solve numerically the equations

$\displaystyle\dot{t}=\frac{E}{\Lambda^{2}\,\left(1-\frac{2M}{r}\right)},$
(67) $\displaystyle\dot{\phi}=\frac{L\,\Lambda^{2}}{r^{2}\sin^{2}\theta},$
(68)

obtained from Eq. (20), together with

$\displaystyle\ddot{r}+\Gamma^{r}_{\ \mu\nu}\dot{x}^{\mu}\dot{x}^{\nu}=0,$
(69) $\displaystyle\ddot{\theta}+\Gamma^{\theta}_{\
\mu\nu}\dot{x}^{\mu}\dot{x}^{\nu}=0,$ (70)

obtained from the geodesic equation. We need to supply initial conditions for
Eqs. (67)-(70). This is achieved by computing the 4-momentum of the photon, as
measured by a static observer. The tetrad basis of this static observer in the
SMBH spacetime can be written as

$\displaystyle\hat{e}^{\hat{0}}_{\
\mu}=\left(\sqrt{1-\frac{2M}{r}}\Lambda,0,0,0\right),$ (71)
$\displaystyle\hat{e}^{\hat{1}}_{\
\mu}=\left(0,\frac{\Lambda}{\sqrt{1-\frac{2M}{r}}},0,0\right),$ (72)
$\displaystyle\hat{e}^{\hat{2}}_{\ \mu}=\left(0,0,r\,\Lambda,0\right),$ (73)
$\displaystyle\hat{e}^{\hat{3}}_{\
\mu}=\left(0,0,0,\frac{r\sin\theta}{\Lambda}\right),$ (74)

obeying

$\hat{e}^{\hat{a}}_{\ \mu}\hat{e}^{\hat{b}}_{\
\nu}n_{\hat{a}\hat{b}}=g_{\mu\nu},$ (75)

where $n_{\hat{a}\hat{b}}$ are the covariant components of the Minkowski
metric tensor. $\hat{e}^{\hat{0}}_{\ \mu}$ is the 4-velocity of the static
observer, and $\hat{e}^{\hat{i}}_{\ \mu}$ (i=1, 2, 3) are the spatial
directions of the reference frame. The components of the 4-momentum of the
photon measured in the static frame, namely

$p^{\hat{a}}=\hat{e}^{\hat{a}}_{\ \mu}\,p^{\mu},$ (76)

are given by

$\displaystyle p^{\hat{t}}=E_{obs}=\frac{E}{\sqrt{1-\frac{2M}{r}}\Lambda},$
(77) $\displaystyle p^{\hat{r}}=\frac{\Lambda}{\sqrt{1-\frac{2M}{r}}}\dot{r},$
(78) $\displaystyle p^{\hat{\theta}}=r\Lambda\dot{\theta},$ (79)
$\displaystyle p^{\hat{\phi}}=\frac{L\Lambda}{r\sin\theta}.$ (80)

The photon’s linear 3-momentum $\vec{p}$ in the static frame has the
components

$\vec{p}=\left(p^{\hat{r}},p^{\hat{\theta}},p^{\hat{\phi}}\right),$ (81)

which are parametrized as functions of the celestial coordinates
$\left(\alpha,\beta\right)$ as

$\displaystyle p^{\hat{r}}=|\vec{p}|\cos\alpha\cos\beta,$ (82) $\displaystyle
p^{\hat{\theta}}=|\vec{p}|\sin\alpha,$ (83) $\displaystyle
p^{\hat{\phi}}=|\vec{p}|\cos\alpha\sin\beta.$ (84)

From Eqs. (77)-(80) and Eqs. (82)-(84), we obtain

$\displaystyle
E=|\vec{p}|\left.\Lambda\sqrt{1-\frac{2M}{r}}\right|_{(r_{\textrm{obs}},\theta_{\textrm{obs}})},$
(85)
$\displaystyle\dot{r}=|\vec{p}|\Lambda^{-1}\left.\sqrt{1-\frac{2M}{r}}\cos\alpha\cos\beta\right|_{(r_{\textrm{obs}},\theta_{\textrm{obs}})},$
(86)
$\displaystyle\dot{\theta}=|\vec{p}|\left.\Lambda^{-1}\frac{\sin\alpha}{r}\right|_{(r_{\textrm{obs}},\theta_{\textrm{obs}})},$
(87) $\displaystyle
L=|\vec{p}|\left.\Lambda^{-1}r\sin\theta\cos\alpha\sin\beta\right|_{(r_{\textrm{obs}},\theta_{\textrm{obs}})},$
(88)

where $(r_{\textrm{obs}},\theta_{\textrm{obs}})$ are the observer’s
coordinates. The observer’s position matches the initial photons’ position,
and Eqs. (85)-(88) are the initial conditions for the equations of motion
(67)-(70). Each pair $(\alpha,\beta)$ represents one point in the observer’s
local sky. The shadow corresponds to the points in the $\alpha$ – $\beta$
plane for which the light rays, evolved backwards in time, using a Dormand-
Prince method [80], are captured by the BH. We compared our results, obtained
with a C++ code with the ones obtained using the pyhole code [81] and obtained
excellent agreement. In order to study the behavior of the scattered light
rays, we place a celestial “sphere" at $r_{cs}>r_{\textrm{obs}}$, with a
different color in each quadrant (red, green, blue, yellow). We also place a
white circular spot on the celestial sphere along the axis joining the
observer and the origin of the radial coordinate — see Fig. 6 top panels.
Unlike asymptotically flat localized objects, where all interesting lensing
features occur near a gravitational center, in asymptotically Melvin
spacetimes nontrivial features may appear panoramically. Thus, it is of
interest to consider panoramic images, with a 360$\degree$ horizontal and
180$\degree$ vertical span. When the BH and the magnetic field are absent
(i.e. $B=0$, $M=0$), the observer’s panoramic view is shown in Fig. 6 (bottom
panel).

Figure 6: Top left: the observer’s local sky in the equatorial position with
reference directions $\\{N,S,F,B,E,W\\}$. Top right: the light emitting
colored celestial sphere encircling both the observer and the BH. Bottom: the
panoramic image of the observer’s local sky in flat spacetime.

(a) Minkowski spacetime ($B=0=M$), with $r_{p}=1$, $r_{cs}=2.5$ and
$\theta_{\textrm{obs}}=\pi/2$.

(b) Melvin spacetime ($B\neq 0=M$), with $B\,r_{p}=1$, $B\,r_{cs}=2.5$ and
$\theta_{\textrm{obs}}=\pi/2$.

Figure 7: Panoramic image obtained in Minkowski (upper panel) and Melvin
(lower panel) spacetimes. The latter exhibits gravitational lensing.

In the following, SMBHs ($M,B\neq 0$) and Melvin universes ($M=0,B\neq 0$)
will be said to be observed under comparable conditions if the perimetral
radius $r_{p}$, given by

$r_{p}\equiv\left.\sqrt{g_{\phi\phi}}\right|_{\theta=\pi/2}=\frac{4r}{4+B^{2}r^{2}},$
(89)

is the same in both cases (rather than the radial coordinate
$r_{\textrm{obs}}$). In contrast to Schwarzschild [8], if $B\neq 0$ the
perimetral radius has a maximum value, given by

$r_{p}^{max}=\frac{1}{B},$ (90)

and the corresponding value of the radial coordinate is

$r=\frac{2}{B}.$ (91)

Since $r\in\ ]0,+\infty[$, for SMBH spacetimes with $B\neq 0$ there will be
generically two radial coordinates $r$ for the same perimetral radius $r_{p}$.
The observer location is chosen to be (of the two) the one with smallest $r$
value.

###  IV.2 Gravitational lensing in the Melvin universe

The gravitational lensing in the Melvin universe ($B\neq 0=M$) is exhibited in
Fig. 7(b), where it is compared with the Minkowski case, Fig. 7(a). For the
Melvin universe, $B$ is the only scale and so all dimensionful quantities are
normalized by $B$. The panoramic image is for an observer at $r_{p}=1/B$,
$\theta_{\textrm{obs}}=\pi/2$ and a celestial sphere at $r_{cs}=2.5/B$.

The Melvin spacetime (unlike Minkowski), exhibits gravitational lensing
leading to several notable features:

1)

Looking forward (toward F), the observer sees the celestial sphere flipped
horizontally, i.e., left-right inverted. This left-right flipping is a
consequence of the axial symmetry together with staticity. By contrast, in
spherical symmetry, a gravitational centre tends to flip the image with
respect to the symmetry center, rather than an axis — see, $e.g.,$ Fig. 3 (top
right) in [15].

2)

Looking east/west (toward E/W), the observer sees some chaotic regions. An
understanding of these regions can be obtained by following illustrative
photons. For this purpose we have highlighted two points (1 and 2) of the
image within the turbulent regions. Introducing a radial coordinate that
compactifies the radial direction

$\mathcal{R}_{M}\equiv\frac{r}{1+r}.$ (92)

Figure 8: Left: contour plots of the effective potential for the Melvin
spacetime in the $\mathcal{R}_{M}-\theta$ plane, for the points marked in Fig.
7(b). The trajectories start off at the blue dot (observer) and terminate at
the black dot, where they meet the celestial sphere (red dotted vertical
line). Right: trajectories described by the photon for the corresponding
values of $\eta$ (displayed in the left column).

Figure 8 exhibits the effective potential $H(r,\theta)$ contour plot, in terms
of $\mathcal{R}_{M}$, for the Melvin universe and the photons corresponding to
the points highlighted in Fig. 7(b). The black regions in the left panels of
Fig. 8 are the inaccessible regions for the photon due to the restriction (31)
on the effective potential. Recall that spatial infinity is inaccessible to
these photons (or any photons that have $\eta\neq 0$). The photons start off
at the blue dot (observer) and terminate their trajectories at the black dot,
where they meet the celestial sphere (red dotted line). The left panels of
Fig. 8, clarify that the turbulent lensing region arises due to the existence
of multiple turning points of the orbit before encountering the celestial
sphere, which occurs at larger latitudes than that of the initial observation.

3)

Looking north/south or back (toward N/S or B), the observer sees a
qualitatively similar image as to that in Minkowski spacetime in terms of the
color patterns, albeit a different structure for the grid lines, likely a
consequence of $r$=constant surfaces in (II) not being round spheres.

###  IV.3 Undercritical SMBHs

We now turn our attention to the general SMBH spacetime ($M\neq 0$)
considering not only the lensing but also the shadows produced by the
existence of a BH region. Let us start by considering how the Schwarzschild
shadow is deformed due to small magnetic fields and undercritical SMBHs. For
this purpose, we shall consider nonpanoramic, narrrow angle images, centered
around the F direction.

Figure 9 exhibits the shadows and lensing of SMBHs for $BM\ll 1$, which is
(likely) the parameter range most suitable to model real astrophysical
environments (near the BH region). In this figure the observer is located on
the equatorial plane ($\theta_{\textrm{obs}}=\pi/2$), at $r_{p}=10M$ and the
viewing angle is $75\degree$ (as opposed to $360\degree$ in the panoramic
images). The colored celestial sphere in Fig. 9 has radius $r_{cs}=25M$. For
$B=0$, we recover the Schwarzschild BH shadow, which is perfectly circular in
the observer’s local sky. For $B>0$, we note that the shadow is no longer
circular, becoming oblate (elongated along the horizontal direction) as the
value of $B$ is increased. It is a rather curious feature in the SMBH
spacetime that the intrinsic horizon geometry becomes prolate (cf. Fig.1) but
the BH shadow becomes oblate (cf. Fig. 9). This is yet another manifestation
of how the shadow is not a faithful probe of the event horizon geometry [82].

(a) BM=0

(b) BM=0.03

(c) BM=0.04

(d) BM=0.05

Figure 9: Nonpanoramic images of the shadows and gravitational lensing of
(undercritical) SMBHs for different values of the magnetic field $B$. We have
chosen the observer perimetral radius $r_{p}=10M$,
$\theta_{\textrm{obs}}=\pi/2$ and $r_{cs}=25M$.

Figure 10: Top: panoramic image of the shadow and gravitational lensing of an
overcritical SMBH ($BM=0.2$). The observer is located at $r_{p}=4.698M$,
$\theta_{\textrm{obs}}=\pi/2$ and $r_{cs}=45M$. Middle: the celestial sphere
was colored white to emphasize the multiple disconnected shadows (black
regions) in the top image. Bottom: three FPOs (blue, green and black curves),
corresponding to three points at the shadow edge in the middle panel, at the
intersection of vertical and horizontal lines with the corresponding color.

Further analyzing the impact of B on the shadow in Fig. 9, we remark that in
the Schwarzschild case, the white spot on the celestial sphere appears as a
ring in the observer’s local sky: an Einstein ring [83]. For the SMBH, on the
other hand, the spherically symmetric Einstein ring is broken up; as the value
of $B$ is increased one observes different patterns. For instance, for
$BM=0.04$ [Fig. 9(c)] four (distorted) copies of the white spot are seen,
forming the ends of a cross centered in the BH shadow, somewhat reminiscent of
the famous Einstein cross formed by the lensed image of the quasar QSO
2237+0305.

###  IV.4 Overcritical SMBHs

Let us now consider large magnetic fields and overcritical SMBHs. As we have
seen in Sec. III.1, for $B>B_{c}$ one of the most remarkable features
concerning null geodesics in the SMBH spacetime arises: there are no LRs. This
leads to the following intriguing question: since LRs determine the edge of
the BH shadow (on the equatorial plane), what happens to the shadow of SMBHs
with $B>B_{c}$?

To investigate this question, Fig. 10 exhibits the shadow and gravitational
lensing for a SMBH with $BM=0.2$ [thus $B>B_{c}$, cf. Eq. (45)]. This is a
panoramic image since the impact of the absence of LRs is to permit the
possibility of a panoramic, (almost) $360\degree$ wide, BH shadow.

Let us dissect the previous statement. In Fig. 10, the observer’s perimetral
radius and the radius of the celestial sphere are $r_{p}=4.698M$ and
$r_{cs}=45M$, respectively. The location of the celestial sphere impacts on
how panoramic the shadow is. Figure 2 makes clear that outwards directed
equatorial light rays with $\eta\neq 0$ always reach a radial turning point.
This turning point is at an increasingly large $r$ when approaching the
backwards (B) direction, and $r\rightarrow\infty$ as $\eta=0$, i.e. as the
observer looks toward B. Hitting the celestial sphere is a stoppage criterion,
yielding a colorful pixel. But all lights rays except the one with $\eta=0$,
would reach a radial turning point whence they necessarily fall back into the
BH region (green curve in Fig. 2) if the celestial sphere is pushed to a
sufficiently large $r_{cs}$, yielding a panoramic (almost) $360\degree$ wide
BH shadow.444Strictly speaking, there is always a vicinity of the B direction
(backwards), wherein no shadow is seen, since such photons hit the celestial
sphere before reaching their radial turning point no matter how large (but
finite) $r_{cs}$ is. But for a sufficiently far away celestial sphere this
vicinity is reduced, in the discretized image, to a single point: the B
direction.

Besides the panoramic (equatorial) shadow, the image of overcritical SMBHs in
Fig. 10 has other outstanding features. In Fig. 10 the shadow image is
exhibited with a white colored celestial sphere, to appreciate the appearance
of multiple disconnected shadows. Thus both the shadow and the lensing [as
seen in Fig. 10] contain imprints of a chaotic pattern.

In order to better understand these chaotic patterns in the overcritical SMBH
we proceed as in the Melvin universe case, analysing individual orbits. In
Fig. 10 we have highlighted three distinct points on the celestial sphere: one
scattered orbit in the chaotic region (point 1) and two absorbed orbits in two
distinct shadow regions (points 2 and 3). The analysis of the corresponding
orbits is done in Fig. 11, specifying the value of $\eta$ and the contour
plots of the effective potential $H(r,\theta)$. As in Fig. 8, the blue (black)
dots mark the beginning (end) of the orbit, corresponding to the position of
the observer (celestial sphere or apparent horizon). In the right panels of
Fig. 11, the black sphere has radius equal to the one of the BH apparent
horizon. In the left panels of Fig. 11, the radial coordinate along the
horizontal axis is the compactified coordinate [81]

$\displaystyle\mathcal{R}\equiv\frac{\sqrt{r^{2}-(2M)^{2}}}{1+\sqrt{r^{2}-(2M)^{2}}},$
(93)

such that spatial infinity corresponds to $\mathcal{R}=1$, and the apparent
horizon to $\mathcal{R}=0$.

Figure 11: Left: contour plots of the effective potential in the
$\mathcal{R}-\theta$ plane, for the points selected in Fig. 10. Right:
trajectories described by the photon for the corresponding values of $\eta$
(displayed in the left panels).

The photon trajectory corresponding to point 1 shows how a north directed
observation sees the southern hemisphere, where the photon reaches the
celestial sphere. On the other hand, trajectory 2 shows that a slightly
different observation direction is described by a photon that reaches the
radial turning point (rather than the celestial sphere), whence it falls into
the BH region. The same occurs for point 3 whose effective potential pushes
the whole of the celestial sphere into the forbidden region, explaining why no
colorful points can be found in its vicinity. Besides the $r-\theta$ motion,
Fig. 11 also exhibits the nontrivial azimuthal motion, which adds to the
possible color mixing in the lensing regions with chaotic patterns.

Before closing the analysis of Fig. 10 there is one last point that merits
discussion. The absence of equatorial LRs determines there is no shadow edge
on the equatorial plane, apart from the one determined by the arbitrary
placement of the celestial sphere. How about outside the equatorial plane? In
asymptotically flat cases, the BH shadow edge in nonspherically symmetric
spacetimes is determined, generically, by a set of bound photon orbits, dubbed
fundamental photon orbits (FPOs) in [84]. For the Kerr case these are the
familiar spherical orbits, which degenerate in the corotating and
counterrotating LRs on the equatorial plane. Are there nonequatorial FPOs in
the SMBH spacetime?

The answer to this question is yes. Since the SMBH does not admit
separability, we have numerically computed such FPOs. These are illustrated in
Fig. 10, for the overcritical SMBH with $BM=0.2$. The corresponding points in
the shadow edge, for each of the three illustrated FPOs, are identified by the
intersection between the vertical and horizontal dashed lines in Fig. 10 (with
the same color). It is worth remarking that the black FPO is planar; but it is
polar rather than equatorial, and it is geometrically noncircular, as we now
explain. First, observe that identical polar planar FPOs occur for any
$\phi=\text{constant}$, due to the axial symmetry and always intersect the $z$
axis. Second, inspection of the $\phi=\text{constant}$ hypersurfaces of the
SMBH geometry (II) shows that these sections are conformal to the
Schwarzschild ones:

$\displaystyle\left.ds^{2}\right|_{\phi=\text{cte.}}=\Lambda^{2}\left[-\left(1-\frac{2\,M}{r}\right)dt^{2}+\frac{dr^{2}}{\left(1-\frac{2\,M}{r}\right)}+r^{2}d\theta^{2}\right].$
(94)

Thus the usual Schwarzschild LR at $r=3M$ remains a geodesic here, but now the
$\theta$ motion plays the role of the usual $\phi$ in Schwarzschild and thus
these are polar rather than equatorial orbits. These are our polar planar
FPOs. Third, since these polar orbits with $r=3M$ exist for every $\phi$, they
span the $r=3M$, $t=constant$, 2-surface of the SMBH geometry (II):

$\displaystyle\left.ds^{2}\right|_{t=\text{cte.},\ r=3M}=$
$\displaystyle(3M)^{2}\left[\Lambda^{2}(r,\theta)\,d\theta^{2}+\frac{\,\sin^{2}\theta}{\Lambda^{2}(r,\theta)}\,d\phi^{2}\right].$
(95)

This (topologically) 2-sphere is not round (but rather prolate, as the
apparent horizon, Fig. 1). It has equatorial $\ell_{\rm eq}$ and polar
$\ell_{\rm p}$ proper lengths:

$\displaystyle\frac{\ell_{\rm eq}}{2\pi(3M)}=$
$\displaystyle\frac{1}{1+\frac{(3MB)^{2}}{4}},$ (96)
$\displaystyle\frac{\ell_{\rm p}}{2\pi(3M)}=$ $\displaystyle
1+\frac{(3MB)^{2}}{8}.$ (97)

Thus, we conclude that, in this sense, for $B\neq 0$ the polar FPOs are not
circular, but rather elongated in the magnetic field direction, since
$\ell_{\rm p}>\ell_{\rm eq}$. Similar planar FPOs can be found in examples of
di-BHs — see [85].

The existence of FPOs suggests that there are shadow edges that are
independent of the celestial sphere location (unlike, say, the equatorial
shadow edge). We also remark that FPOs also exists for the undercritical SMBHs
and they determine the BH shadow edges shown in Fig. 9.

##  V Final remarks

The SMBH spacetime first obtained by Ernst in 1976 [53] as a solution of the
Einstein-Maxwell (electrovacuum) equations, is interpreted as a Schwarzschild
BH in a Melvin magnetic Universe. The latter is a sort of confining box, for
light rays, that can only reach infinity if they have no angular momentum. The
analysis of the null geodesics in this paper has shown that the SMBH provides
a rich and novel case study for features of gravitational lensing and BH
shadows. As a summary of our results we would like to emphasize the following
features:

$\bullet$

The empty Melvin universe ($M=0\neq B$) [54] has an equatorial stable LR,
which, moreover, due to a “vertical" isometry, implies the existence of a tube
of planar LRs.

$\bullet$

The SMBH family ($M\neq 0\neq B$) yields a novel example of BH spacetimes
wherein the LRs total topological charge [5; 4] is zero. This can be
interpreted as the cancellation between the $w=-1$ Schwarzschild and the
$w=+1$ Melvin (asymptotics) topological charges.

$\bullet$

The vanishing of the topological LR charge permits that, unlike for
asymptotically flat [4] and asymptotically de Sitter or AdS [6] BHs, there can
be asymptotically Melvin BHs without LRs.

$\bullet$

The SMBHs family is then naturally divided into undercritical ($B<B_{c}$) and
overcritical ($B>B_{c}$) spacetimes, where $B_{c}$ is given by Eq. (45). The
latter describe BH spacetimes (regular on and outside a horizon) without LRs
outside the horizon.

$\bullet$

Performing ray-tracing, we have shown that the absence of equatorial LRs
yields a panoramic shadow: even a far away observer will see an (almost)
$360\degree$ BH shadow, looking at any direction along the equatorial plane,
as long as the celestial sphere is sufficiently far away.

$\bullet$

Despite the absence of LRs, there are nonequatorial (unstable) FPOs, that
determine parts of the shadow edge, outside the equatorial plane. There are,
in particular planar, polar (but noncircular) photon orbits, which determine
the north and south pole of the shadow edge.

$\bullet$

Despite their simplicity, both the Melvin universe and the MSBH family yield
chaotic lensing (and the latter also multiple shadows). Some insight into this
behavior can be obtained by following illustrative orbits and analyzing the
corresponding effective potentials.

$\bullet$

For small $BM$ (undercritical) SMBHs, the shadow becomes oblate, in a curious
contrast with the intrinsic horizon geometry which becomes prolate.

Although most of our results are unlikely to be astrophysically relevant
(since Melvin asymptotics do not describe the real Universe), the shadow
deformation for small magnetic fields could, potentially, describe those for
an astrophysical BH immersed in a sufficiently strong poloidal magnetic field
whose backreaction would have to be taken into account. Interestingly, the
recent EHT results suggest that the M87* magnetic field is poloidal, rather
than toroidal [52].

Finally, our work suggests further studies. For instance, a deeper
understanding of the chaotic patterns (and multiple shadows) as well as their
dependence on the celestial sphere would be desirable. These features are
related to the existence of closed pockets and forbidden regions in the
effective potential, leading to multiple turning points in the photons’
trajectories. We have observed, for instance, that the time that a photon
takes to reach the celestial sphere increases in the chaotic region due to the
existence of several radial turning points, which has also been observed in
other examples [15; 81; 85; 86]. As a second possible extension, a natural
follow up with be to consider either the Reissner-Nordström-Melvin or the
Kerr-Melvin solutions (see, $e.g.,$ [87]). Both of these acquire angular
momentum, which will lead to frame-dragging effects in the lensing and
shadows.

###### Acknowledgements.

The authors thank Fundação Amazônia de Amparo a Estudos e Pesquisas (FAPESPA),
Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq) and
Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (Capes) - Finance
Code 001, in Brazil, for partial financial support. This work is supported by
the Center for Research and Development in Mathematics and Applications
(CIDMA) through the Portuguese Foundation for Science and Technology (FCT -
Fundação para a Ciência e a Tecnologia), references No. UIDB/04106/2020, No.
UIDP/04106/2020, and No. BIPD/UI97/7484/2020. We acknowledge support from the
projects No. PTDC/FIS-OUT/28407/2017, No. CERN/FIS-PAR/0027/2019, and No.
PTDC/FIS-AST/3041/2020. This work has further been supported by the European
Union’s Horizon 2020 research and innovation (RISE) programme H2020-MSCA-
RISE-2017 Grant No. FunFiCO-777740. The authors would like to acknowledge
networking support by the COST Action CA16104.

## References

  * [1] L. C. B. Crispino and S. Paolantonio, The first attempts to measure light deflection by the Sun, Nat. Astron. 4, 6 (2020).
  * [2] F. W. Dyson, A. S. Eddington and C. Davidson, A determination of the deflection of light by the Sun’s gravitational field, from observations made at the total eclipse of May 29, 1919, Philos. Trans. R. Soc. Lond. A 220, 291 (1920).
  * [3] L. C. B. Crispino and D. J. Kennefick, A hundred years of the first experimental test of general relativity, Nat. Phys. 15, 416 (2019).
  * [4] P. V. P. Cunha and C. A. R. Herdeiro, Stationary Black Holes and Light Rings, Phys. Rev. Lett. 124, 181101 (2020).
  * [5] P. V. P. Cunha, E. Berti and C. A. R. Herdeiro, Light-Ring Stability for Ultracompact Objects, Phys. Rev. Lett. 119, 251102 (2017).
  * [6] S. W. Wei, Topological charge and black hole photon spheres, Phys. Rev. D 102, 064039 (2020).
  * [7] M. Guo and S. Gao, Universal properties of light rings for stationary axisymmetric spacetimes, Phys. Rev. D 103, 104031 (2021).
  * [8] P. V. P. Cunha and C. A. R. Herdeiro, Shadows and strong gravitational lensing: A brief review, Gen. Relativ. Gravit., 50, 42 (2018).
  * [9] H. Falcke, F. Melia, and E. Agol, Viewing the shadow of the black hole at the galactic center, Astrophys. J. Lett. 528, L13 (2000)
  * [10] R. P. Kerr, Gravitational Field of a Spinning Mass as an Example of Algebraically Special Metrics, Phys. Rev. Lett. 11, 237(1963).
  * [11] J. M. Bardeen, Timelike and null geodesics in the Kerr metric, in Black Holes (Les Astres Occlus) (edited by C. Dewitt and B. S. Dewitt), (1973), pp. 215-219.
  * [12] The Event Horizon Telescope Collaboration, First M87 event horizon telescope results. I. The shadow of the supermassive black hole, Astrophys, J. Lett. 875, L1 (2019).
  * [13] A. de Vries, The apparent shape of a rotating charged black hole, closed photon orbits and the bifurcation set $A_{4}$, Classical Quantum Gravity 17, 123 (2000).
  * [14] N. Tsukamoto, Black hole shadow in an asymptotically flat, stationary, and axisymmetric spacetime: The Kerr-Newman and rotating regular black holes, Phys. Rev. D 97, 064021 (2018).
  * [15] P. V. P. Cunha, C. A. R. Herdeiro, E. Radu, and H. F. Rúnarsson, Shadows of Kerr Black Holes with Scalar Hair, Phys. Rev. Lett. 115, 211102 (2015).
  * [16] P. V. P. Cunha, C. A. R. Herdeiro, and E. Radu, EHT constraint on the ultralight scalar hair of the M87 supermassive black hole, Universe 5, 220 (2019).
  * [17] P. V. P. Cunha, J. Grover, C. Herdeiro, E. Radu, H. Runarsson, and A. Wittig, Chaotic lensing around boson stars and Kerr black holes with scalar hair, Phys. Rev. D 94, 104023 (2016).
  * [18] K. Hioki and U. Miyamoto, Hidden symmetries, null geodesics, and photon capture in the Sen black hole, Phys. Rev. D 78, 044007 (2008).
  * [19] P. V. P. Cunha, C. A. R. Herdeiro and E. Radu, Spontaneously Scalarized Kerr Black Holes in Extended Scalar-Tensor–Gauss-Bonnet Gravity, Phys. Rev. Lett. 123, 011101 (2019).
  * [20] P. V. P. Cunha, C. A. R. Herdeiro, and E. Radu, Isolated black holes without $\mathbb{Z}_{2}$ isometry, Phys. Rev. D 98, 104060 (2018).
  * [21] P. V. P. Cunha, C. A. R. Herdeiro, and M. J. Rodriguez, Shadows of exact binary black holes, Phys. Rev. D 98, 044053 (2018).
  * [22] P. V. P. Cunha, C. A. R. Herdeiro, B. Kleihaus, J. Kunz, and E. Radu, Shadows of Einstein–dilaton–Gauss–Bonnet black holes, Phys. Lett. B 768, 373 (2017).
  * [23] F. Bacchini, D. R. Mayerson, B. Ripperda, J. Davelaar, H. Olivares, T. Hertog, and B. Vercnocke, Fuzzball shadows: Emergent horizons from microstructure, arXiv:2103.12075.
  * [24] A. Eichhorn and A. Held, From a locality-principle for new physics to image features of regular spinning black holes with disks, J. Cosmol. Astropart. Phys. 05 (2021) 073.
  * [25] M. Afrin, R. Kumar, and S. G. Ghosh, Parameter estimation of hairy Kerr black holes from its shadow and constraints from M87, Mon. Not. R. Astron. Soc. 504, 5927 (2021).
  * [26] K. Glampedakis and G. Pappas, Can supermassive black hole shadows test the Kerr metric?, arXiv:2102.13573.
  * [27] M. Zhang and J. Jiang, NUT charges and black hole shadows, Phys. Lett. B 816, 136213 (2021).
  * [28] H. C. D. L. Junior, L. C. B. Crispino, P. V. P. Cunha, and C. A. R. Herdeiro, Can different black holes cast the same shadow?, Phys. Rev. D 103, 084040 (2021).
  * [29] W. L. Qian, S. Chen, C. G. Shao, B. Wang, and R. H. Yue, Cuspy and fractured black hole shadows in a toy model with axisymmetry, arXiv:2102.03820.
  * [30] C. A. R. Herdeiro, A. M. Pombo, E. Radu, P. V. P. Cunha, and N. Sanchis-Gual, The imitation game: Proca stars that can mimic the Schwarzschild shadow, J. Cosmol. Astropart. Phys. 04, (2021) 051.
  * [31] B. H. Lee, W. Lee, and Y. S. Myung, Shadow cast by a rotating black hole with anisotropic matter, Phys. Rev. D 103, 064026 (2021).
  * [32] S. H. Völkel, E. Barausse, N. Franchini, and A. E. Broderick, EHT tests of the strong-field regime of General Relativity, arXiv:2011.06812.
  * [33] M. Zhang and J. Jiang, Shadows of accelerating black holes, Phys. Rev. D 103, 025005 (2021).
  * [34] E. Contreras, Á. Rincón, G. Panotopoulos, and P. Bargueño, Geodesic analysis and black hole shadows on a general non–extremal rotating black hole in five–dimensional gauged supergravity, arXiv:2010.03734.
  * [35] T. C. Ma, H. X. Zhang, H. R. Zhang, Y. Chen, and J. B. Deng, Shadow cast by a rotating and nonlinear magnetic-charged black hole in perfect fluid dark matter, Mod. Phys. Lett. A 36, 2150112 (2021).
  * [36] F. Long, S. Chen, M. Wang, and J. Jing, Shadow of a disformal Kerr black hole in quadratic degenerate higher-order scalar–tensor theories, Eur. Phys. J. C 80, 1180 (2020).
  * [37] A. Belhaj, M. Benali, A. E. Balali, W. E. Hadri, H. El Moumni, and E. Torrente-Lujan, Black hole shadows in M-theory scenarios, Int. J. Mod. Phys. D 30, 2150026 (2021).
  * [38] V. I. Dokuchaev and N. O. Nazarova, Visible shapes of black holes M87* and SgrA*, Universe 6, 154 (2020).
  * [39] M. Khodadi, A. Allahyari, S. Vagnozzi, and D. F. Mota, Black holes with scalar hair in light of the Event Horizon Telescope, J. Cosmol. Astropart. Phys 09, 026 (2020).
  * [40] J. Badía and E. F. Eiroa, Influence of an anisotropic matter field on the shadow of a rotating black hole, Phys. Rev. D 102, 024066 (2020).
  * [41] A. Övgün and İ. Sakallı, Testing generalized Einstein–Cartan–Kibble–Sciama gravity using weak deflection angle and shadow cast, Classical Quantum Gravity 37, 225003 (2020).
  * [42] S. Chen, M. Wang, and J. Jing, Polarization effects in Kerr black hole shadow due to the coupling between photon and bumblebee field, J. High Energy Phys. 07 (2020) 054.
  * [43] G. Creci, S. Vandoren, and H. Witek, Evolution of black hole shadows from superradiance, Phys. Rev. D 101, 124051 (2020).
  * [44] A. Abdujabbarov, M. Amir, B. Ahmedov, and S. G. Ghosh, Shadow of rotating regular black holes, Phys. Rev. D 93, 104004 (2016).
  * [45] S. W. Wei and Y. X. Liu, Testing the nature of Gauss-Bonnet gravity by four-dimensional rotating black hole shadow, Eur. Phys. J. Plus 136, 436 (2021).
  * [46] F. H. Vincent, Z. Meliani, P. Grandclement, E. Gourgoulhon, and O. Straub, Imaging a boson star at the Galactic center, Classical Quantum Gravity 33, 105015 (2016).
  * [47] R. D. Blandford and R. L. Znajek, Electromagnetic extraction of energy from Kerr black holes. Mon. Not. R. Astron. Soc. 179, 433 (1977).
  * [48] K. Mori et al., NuSTAR discovery of a 3.76s transient magnetar near Sagittarius A*. Astrophys. J. Lett. 770, L23 (2013).
  * [49] J. A. Kennea et al., SWIFT discovery of a new soft gamma repeater, SGR J1745-29, near Sagittarius A*, Astrophys. J. Lett. 770, L24 (2013).
  * [50] R. P. Eatough et al., A strong magnetic field around the supermassive black hole at the centre of the Galaxy, Nature (London) 501, 391 (2013).
  * [51] S. A. Olausen and V. M. Kaspi, The McGill magnetar catalog. Astrophys. J. Suppl. Ser. 212, 6 (2014).
  * [52] K. Akiyama et al. (The Event Horizon Telescope Collaboration), First M87 event horizon telescope results. VIII. Magnetic field structure near the event horizon, Astrophys. J. Lett. 910, L13 (2021).
  * [53] F. J. Ernst, Black holes in a magnetic universe. J. Math. Phys. (N. Y.) 17, 54 (1976).
  * [54] M. A. Melvin, Pure magnetic and electric geons, Phys. Lett. 8, 65 (1964).
  * [55] N. Dadhich, C. Hoenselaers, and C. V. Vishveshwara, Trajectories of charged particles in the static Ernst space-time, J. Phys. A 12 215 (1979).
  * [56] V. Karas and D. Vokrouhlický, Chaotic motion of test particles in the Ernst space-time, Gen. Relativ. Gravit. 24, 7 (1992).
  * [57] D. Vokrouhlický and V. Karas, Motion of charged particles around a black hole in a magnetic field, Astron. Astrophys. 243, 165 (1991).
  * [58] Y. K. Lim, Motion of charged particles around a magnetized/electrified black hole, Phys. Rev. D 91, 024048 (2015).
  * [59] E. P. Esteban, Geodesics in the Ernst metric, Il Nuovo Cimento 79, 1 (1984).
  * [60] S. V. Dhurandhar and D. N. Sharma, Null geodesics in the static Ernst space-time, J. Phys. A 16, 99 (1983).
  * [61] Z. Stuchlík and S. Hledík, Photon capture cones and embedding diagrams of the Ernst spacetime. Classical Quantum Gravity 16, 1377 (1999).
  * [62] R. M. Wald, Black hole in a uniform magnetic field, Phys. Rev. D 10, 1680 (1974).
  * [63] J. Chen, H. Liao, Y. Wang, and T. Chen, Absorption and scattering of scalar wave from Schwarzschild black hole surrounded by magnetic field, Eur. Phys. J. C 73, 2395 (2013).
  * [64] R. A. Konoplya and R. D. B. Fontana, Quasinormal modes of black holes immersed in a strong magnetic field, Phys. Lett. B 659, 375 (2008).
  * [65] R. Brito, V. Cardoso, and P. Pani, Superradiant instability of black holes immersed in a magnetic field, Phys. Rev. D 89, 104045 (2014).
  * [66] R. A. Konoplya, Superradiant instability for black holes immersed in a magnetic field, Phys. Lett. B 666, 283 (2008); 670, 459 (E) (2009).
  * [67] N. M. Santos and C. A. R. Herdeiro, Black holes, stationary clouds and magnetic fields, Phys. Lett. B 815, 136142 (2021).
  * [68] J. A. Wheeler, Geons, Phys. Rev. 97, 511 (1955).
  * [69] M. A. Melvin and J. S. Wallingford, Orbits in a magnetic universe, J. Math. Phys. (N.Y.) 7, 333 (1966).
  * [70] K. S. Thorne, Absolute stability of Melvin’s magnetic universe, Phys. Rev. 139, B244 (1965).
  * [71] Y. K. Lim, Hypocycloid motion in the Melvin magnetic universe, Phys. Rev. D 101, 104031 (2020).
  * [72] G. W. Gibbons and C. A. R. Herdeiro, The Melvin universe in Born-Infeld theory and other theories of nonlinear electrodynamics, Classical Quantum Gravity 18, 1677 (2001).
  * [73] S. W. Hawking and G. F. R. Ellis, The Large Scale Structure of Space-Time (Cambridge University Press, Cambridge, England, 1973).
  * [74] R. M. Wald, General Relativity (The University of Chicago Press, Chicago, 1984).
  * [75] E. Poisson, A Relativist’s Toolkit: The Mathematics of Black-Hole Mechanics (Cambridge University Press, Cambridge, England, 2004).
  * [76] W. J. Wild and R. M. Kerns, Surface geometry of black hole in a magnetic field, Phys. Rev. D 21, 332 (1980).
  * [77] R. Kulkarni and N. Dadhich, Surface geometry of a rotating black hole in a magnetic field, Phys. Rev. D33, 2780 (1986).
  * [78] G. W. Gibbons, C. A. R. Herdeiro, and C. Rebelo, Global embedding of the Kerr black hole event horizon into hyperbolic 3-space, Phys. Rev. D 80, 044014 (2009).
  * [79] B. Carter, Global structure of the Kerr family of gravitational fields, Phys. Rev. 174, 1559 (1968).
  * [80] W. H. Press, S. A. Teukolsky, and W.T. Vetterling. Numerical Recipes in C++: The Art of Scientific Computing (Cambridge University Press, Cambridge, England, 2002).
  * [81] P. V. P. Cunha, J. Grover, C. Herdeiro, E. Radu, H. Runarsson, and A. Wittig, Chaotic lensing around boson stars and Kerr black holes with scalar hair, Phys. Rev. D 94, 104023 (2016).
  * [82] P. V. P. Cunha, C. A. R. Herdeiro, and M. J. Rodriguez, Does the black hole shadow probe the event horizon geometry?, Phys. Rev. D 97, 084020 (2018).
  * [83] A. Einstein, Science 84, 506 (1936).
  * [84] P. V. P. Cunha, C. A. R. Herdeiro, and E. Radu, Fundamental photon orbits: Black hole shadows and spacetime instabilities, Phys. Rev. D 96, 024039 (2017).
  * [85] J. O. Shipley and S. R. Dolan, Binary black hole shadows, chaotic scattering and the Cantor set, Classical Quantum Gravity 33, 175001 (2016).
  * [86] A. Bohn, W. Throwe, F. Hbert, K. Henriksson, and D. Bunandar, What does a binary black hole merger look like? Classical Quantum Gravity 32, 065002 (2015).
  * [87] G. W. Gibbons, A. H. Mujtaba, and C. N. Pope, Ergoregions in magnetised black hole spacetimes, Classical Quantum Gravity 30, 125008 (2013).


$\Pr\left\\{d(x||U)\leq\epsilon\right\\}\geq\frac{\epsilon}{2}$.

###### Proof.

By Theorem 1 in [6], $d(x||u)\leq\frac{x^{2}}{u}+\frac{(1-x)^{2}}{1-u}-1$.
Therefore, if $u$ satisfies

$u\geq\frac{1}{1+\epsilon}x\quad\text{and}\quad
1-u\geq\frac{1}{1+\epsilon}(1-x)\,,$

then $d(x||u)\leq(1+\epsilon)x+(1+\epsilon)(1-x)-1=\epsilon$. It follows that

$\displaystyle\Pr\left\\{d(x||U)\leq\epsilon\right\\}\geq\Pr\left\\{\frac{1}{1+\epsilon}x\leq
U\leq
1-\frac{1}{1+\epsilon}(1-x)\right\\}=1-\frac{1-x}{1+\epsilon}-\frac{x}{1+\epsilon}=\frac{\epsilon}{1+\epsilon}\geq\frac{\epsilon}{2}\,.$

∎

###### Proof of Proposition 17.

Consider a fixed large $N$. Let $\delta(N)=2N^{-\frac{1}{3}}$. Without loss of
generality, suppose $N$ is large enough such that $\delta(N)<\delta_{0}$ as in
Lemma 18. Let $L(N),s_{0}(N),r_{0}(N)$ be defined for $\delta(N)$ as
$L,s_{0},r_{0}$ are defined for $\delta$ in Lemma 18. If a particle $i$
satisfies that $D_{i}$ is entirely below the line $L(N)$, we say that particle
$i$ is good. Let $E$ be the event that there exists at least one good particle
in $\mathcal{P}_{N}$. It follows that

$\displaystyle\Pr\\{E\\}$ $\displaystyle=1-(1-\Pr\\{\text{particle 1 is
good}\\})^{N}$
$\displaystyle=1-\left(1-\Pr\left\\{d(\theta^{*}_{1}||\theta^{(1)}_{1})\leq\delta(N)\right\\}\cdot\Pr\left\\{d(\theta^{*}_{2}||\theta^{(1)}_{2})\leq\delta(N)\right\\}\right)^{N}$
$\displaystyle\overset{(i)}{\geq}1-\left(1-N^{-1/3}N^{-1/3}\right)^{N}$
$\displaystyle\geq 1-e^{-N^{-2/3}N}=1-e^{-N^{\frac{1}{3}}}\,,$

where $(i)$ is due to
$\Pr\left\\{d(\theta^{*}_{i}||\theta^{(1)}_{i})\leq\delta(N)\right\\}\geq
N^{-\frac{1}{3}}$ by Lemma 19 for $i=1,2$.

Suppose event $E$ is true. Let $i_{0}$ be one good particle. Then by Lemma 18
part (b), any $r\in\mathcal{R}$ must satisfy $r\leq s_{0}(N)$ or $r\geq
r_{0}(N)$. Simple geometry shows that
$s_{0}(N)=\frac{\delta(N)}{\bar{d}_{1}}=\frac{2}{\bar{d}_{1}}N^{-\frac{1}{3}}$
and
$r_{0}(N)=1-\frac{\delta(N)}{\bar{d}_{2}}=1-\frac{2}{\bar{d}_{2}}N^{-\frac{1}{3}}$.
Let $C_{1}=\frac{2}{\bar{d}_{1}}$ and $C_{2}=\frac{2}{\bar{d}_{2}}$, part (a)
of Proposition 17 is proved.

Consider any $\xi\in(s_{0},r_{0})$, let the corresponding dominant particle be
$j$. Then $D_{j}(\xi)\leq D_{i_{0}}(\xi)$. By Lemma 18 part (a), particle $j$
must be action-optimal. Part (b) of Proposition 17 is proved. ∎

###  B.5 Summary

In this section we analyzed $\mathrm{PTS}$ for the two-arm Bernoulli bandit
problem. Our key findings are the following.

  * •

_Fit particles survive, unfit particles decay_ , in the sense described in
Proposition 8 and Conjecture 14. The fitness of a particle $i$ is measured in
terms of its closeness to $\theta^{*}$ by the divergence $D_{i}(r_{t})$, a
convex combination of the KL divergences on the two arms. Unfortunately we
cannot directly compare the fitness of particles because $D_{i}(r_{t})$
depends on the random process $r_{t}$. It is possible that the weights of the
surviving particles oscillates forever due to the counter-reinforcing effect.
Also, the weights of the decaying particles decay exponentially fast.

  * •

_The set of surviving particles is random._ This is mainly due to the self-
reinforcing effect. One way to find out the possible sets of surviving
particles is by drawing the divergence diagram described in Section B.3.

  * •

_Most particles decay_. Under Assumption 3, we expect that all except at most
two particles decay eventually.

  * •

Roughly speaking, with randomly generated particles, _PTS is consistent or
near-consistent with high probability_. See Proposition 16 and Proposition 17.

We believe these findings and some related concepts can be extended to other
and more general kinds of stochastic bandit problems. For example, for the
$K$-arm Bernoulli bandit problem with $K\geq 3$, we expect to observe counter-
reinforcing sets (not just pairs) of particles in PTS, in which the particles
reinforce each other in some way. Proposition 1 provides a generalized method
to identify surviving particles, including counter-reinforcing particles, for
general stochastic bandit problems and for any finite number of particles.

###  B.6 Useful Drift Implied Bounds

This section includes for reference two useful drft implied bounds.

####  B.6.1 One drift implied bound with stochastic dominance

The following result (Proposition 20) is taken out from [9] for convenience of
reference. Let $X_{0},X_{1},\cdots$ be a sequence of random variables. The
drift at time $t$ is defined as $\mathbb{E}[X_{t+1}-X_{t}|\mathcal{F}_{t}]$,
where $\mathcal{F}_{t}=\sigma(X_{0},\cdots,X_{t})$. Consider the following two
conditions:

Condition C1:

$\mathbb{E}\left[(X_{t+1}-X_{t})\mathds{1}_{\\{X_{t}\geq
a\\}}|\mathcal{F}_{t}\right]\leq-\epsilon_{0}\quad t\geq 0$ (22)

for some constants $-\infty\leq a<\infty$ and $\epsilon_{0}>0$. That is, the
drift at time $t$ is strictly negative whenever $X_{t}\geq a$.

Condition C2: There exists a random variable $Z$ with $\mathbb{E}[e^{\lambda
Z}]=D$ for some constants $\lambda>0$ and $D>0$ such that $(\left\lvert
X_{t+1}-X_{t}\right\rvert|\mathcal{F}_{t})\prec Z$. That is, given
$\mathcal{F}_{t}$, $\left\lvert X_{t+1}-X_{t}\right\rvert$ is stochastically
dominated by a random variable with exponential tail.

Let $c,\eta,\rho$ be constants such that

$\displaystyle c$ $\displaystyle\geq\frac{\mathbb{E}[e^{\lambda
Z}]-(1+\lambda\mathbb{E}[Z])}{\lambda^{2}}\,,$ $\displaystyle 0$
$\displaystyle<\eta\leq\lambda\,,$ $\displaystyle\eta$
$\displaystyle<\epsilon_{0}/c\,,$ $\displaystyle\rho$
$\displaystyle=1-\epsilon_{0}\eta+c\eta^{2}\,.$

Then $\rho<1$.

######  Proposition 20 (Theorem 2.3 in [9]).

Conditions C1 and C2 imply that

$P\left\\{X_{t}\geq
b|X_{0}\right\\}\leq\rho^{t}e^{\eta(Y_{0}-b)}+\frac{1-\rho^{t}}{1-\rho}De^{-\eta(b-a)}\,.$

In particular, if $X_{0}\leq a$, then

$P\left\\{X_{t}\geq b|X_{0}\right\\}\leq\frac{D}{1-\rho}e^{-\eta(b-a)}\,.$

####  B.6.2 Another drift implied bound with bounded steps

Two lemmas are stated first.

######  Lemma 21 (Hoeffding’s Lemma).

Suppose $Y$ is a random variable such that $\Pr\left\\{Y\in[a,b]\right\\}=1$,
then
$\mathbb{E}\left[e^{\theta(Y-\mathbb{E}[Y])}\right]\leq\frac{\theta^{2}(b-a)^{2}}{8}$.

######  Lemma 22.

Suppose $(M_{k}:k\geq 0)$ is a non-negative supermartingale. Then for any
$n\geq 0$ and $\gamma>0$, $\Pr\left\\{\max_{0\leq k\leq
n}M_{k}>\gamma\right\\}\leq\frac{\mathbb{E}[M_{0}]}{\gamma}$.

A proof of Lemma 22 can be found in Section 3.4 (Page 69) of [10].

######  Proposition 23.

Consider a randon sequence $(U_{n}:n\geq 1)$ and define
$\mathcal{F}=\varnothing$ and $\mathcal{F}_{k}=\sigma(U_{1},\cdots,U_{k})$.
Suppose $\mathbb{E}[U_{k+1}|\mathcal{F}_{k}]\leq-\mu<0$ for $k\geq 0$ and
$\Pr\left\\{\left\lvert U_{k}\right\rvert\leq C\right\\}=1$ for $k\geq 1$ for
some constancts $\mu,C>0$. Let $X_{n}\triangleq U_{1}+\cdots+U_{n}$ for $n\geq
1$ and $X_{0}=0$. Let $G_{n}\triangleq\max_{0\leq k\leq n}X_{k}$ and
$G\triangleq\max_{k\geq 0}X_{k}$. Then for any $b>0$,
$\Pr\left\\{G>b\right\\}\leq e^{-\frac{2\mu b}{C^{2}}}$ .

###### Proof.

By Hoeffding’s lemma (Lemma 21),

$\mathbb{E}\left[e^{\theta(U_{k}-\mathbb{E}[U_{k}|\mathcal{F}_{k-1}])}|\mathcal{F}_{k-1}\right]\leq
e^{\frac{\theta^{2}(2C)^{2}}{8}}=e^{\frac{\theta^{2}C^{2}}{2}}\,.$

Therefore, for all $k\geq 1$,

$\mathbb{E}\left[e^{\theta U_{k}}|\mathcal{F}_{k-1}\right]\leq
e^{\theta\mathbb{E}\left[U_{k}|\mathcal{F}_{k-1}\right]}e^{\frac{\theta^{2}C^{2}}{2}}\leq
e^{-\theta\mu+\frac{\theta^{2}C^{2}}{2}}\,.$

$-\theta\mu+\theta^{2}C^{2}/2$ is quadratic in $\theta$ and is less than or
equal to zero for all $\theta\in[0,2\mu/C^{2}]$. Let $\theta^{*}=2\mu/C^{2}$.
Then $\mathbb{E}\left[e^{\theta^{*}U_{k}}|\mathcal{F}_{k-1}\right]\leq 1$ for
all $k\geq 1$. Next, define $M_{0}=1$ and $M_{k}=e^{\theta^{*}X_{k}}$ for
$k\geq 1$. $(M_{k}:k\geq 0)$ is a supermartingale because

$\displaystyle\mathbb{E}[M_{k+1}|\mathcal{F}_{k}]=\mathbb{E}\left[e^{\theta^{*}(U_{1}+\cdots+U_{k+1})}|\mathcal{F}_{k}\right]=e^{\theta^{*}(U_{1}+\cdots+U_{k})}\mathbb{E}\left[e^{\theta^{*}U_{k+1}}|\mathcal{F}_{k}\right]=M_{k}\mathbb{E}\left[e^{\theta^{*}U_{k}}|\mathcal{F}_{k-1}\right]\leq
M_{k}\,.$

It follows that, for any $n\geq 0$ and $b>0$,

$\displaystyle\Pr\left\\{G_{n}>b\right\\}$
$\displaystyle=\Pr\left\\{\max_{0\leq k\leq
n}X_{k}>b\right\\}=\Pr\left\\{\max_{0\leq k\leq
n}e^{\theta^{*}X_{k}}>e^{\theta^{*}b}\right\\}=\Pr\left\\{\max_{0\leq k\leq
n}M_{k}>e^{\theta^{*}b}\right\\}$
$\displaystyle\overset{(i)}{\leq}\frac{\mathbb{E}[M_{0}]}{e^{\theta^{*}b}}=e^{-\theta^{*}b}\,.$

Step $(i)$ is due to Lemma 22. Finally, since $G_{n}$ is non-decreasing in $n$
and $G_{n}\rightarrow G$ for each sample path, $\mathds{1}_{\\{G_{n}>b\\}}$ is
non-negative and is non-decreasing in $n$ and
$\mathds{1}_{\\{G_{n}>b\\}}\rightarrow\mathds{1}_{\\{G>b\\}}$ for each sample
path. So by the monotone convergence theorem

$\displaystyle\Pr\left\\{G>b\right\\}=\mathbb{E}\left[\mathds{1}_{\\{G>b\\}}\right]=\lim_{n\rightarrow\infty}\mathbb{E}\left[\mathds{1}_{\\{G_{n}>b\\}}\right]=\lim_{n\rightarrow\infty}\Pr\left\\{G_{n}>b\right\\}\leq
e^{-\theta^{*}b}=e^{-\frac{2\mu b}{C^{2}}}\,.$

∎

######  Corollary 24.

Consider a randon sequence $(U_{n}:n\geq 1)$ and define
$\mathcal{F}=\varnothing$ and $\mathcal{F}_{k}=\sigma(U_{1},\cdots,U_{k})$.
Suppose $\mathbb{E}[U_{k+1}|\mathcal{F}_{k}]\geq\mu>0$ for $k\geq 0$ and
$\Pr\left\\{\left\lvert U_{k}\right\rvert\leq C\right\\}=1$ for $k\geq 1$ for
some constancts $\mu,C>0$. Let $X_{n}\triangleq U_{1}+\cdots+U_{n}$ for $n\geq
1$ and $X_{0}=0$. Let $G_{n}\triangleq\min_{0\leq k\leq n}X_{k}$ and
$G\triangleq\min_{k\geq 0}X_{k}$. Then for any $b>0$, $\Pr\\{G<-b\\}\leq
e^{-\frac{2\mu b}{C^{2}}}$.

###### Proof.

Apply Proposition 23 to the sequence $\\{-X_{n}\\}$. ∎

##  Appendix C Regenerative particle Thompson sampling: choice of hyper-
parameters and more simulations

(a) Bernoulli bandit, $K=10$
$\theta^{*}=[0.05,0.10,\cdots,0.50].$ (b) Bernoulli bandit, $K=100$
$\theta^{*}$ consists of $N=100$ points uniformly spaced over [0.3,0.8]. (c)
Max-Bernoulli bandit, $K=10$, $M=3$
$\theta^{*}=[0.51,0.52,\cdots,0.60].$ (d) Max-Bernoulli bandit, $K=10$, $M=3$
$\theta^{*}=[0.05,0.10,\cdots,0.50].$ (e) Linear bandit, $K=10$,
$\sigma_{W}^{2}=0.1$, $\theta^{*}=[0.2,\cdots,0.2]$. (f) Linear bandit,
$K=100$, $\sigma_{W}^{2}=0.1$, $\theta^{*}=[0.08,\cdots,0.08]$.

Figure 16: More simulations.

The recommended numerical values of the three hyper-parameters for RPTS
(Algorithm 3) are $f_{del}=0.8$, $w_{inact}=0.001$, and $w_{new}=0.01$. The
behavior of the algorithm is relatively insensitive to these values, but
further tuning may be beneficial in a given application. In this section we
comment on how these values influence the performance of the algorithm.

  * •

Analysis for Bernoulli bandits (Section B) and empirical evidence for other
bandit models indicate that with high probability all but a few particles
eventually decay in PTS. Hence it may be attempting to make $f_{del}$ very
large. However, since the set of decaying particles is random, it may happen
that some fit particles end up decaying. Also, a not-so-bad particle may have
an oscillating weight due to counter-reinforcing effects and thus may have low
weight at times. Making $f_{del}$ not too large gives those unfortunate fit
and not-so-bad particles a chance to survive. We have tried $f_{del}=0.8$ and
$f_{del}=0.5$ and both work fine.

  * •

The value of $w_{inact}$ should be small, but if it is too small, it may take
a long time for the CONDITION in Step 9 to become true, especially when the
particles become concentrated in a small subset of the parameter space.

  * •

The value of $w_{new}$ should be small, but strictly larger than $w_{inact}$.
There are three aspects of consideration here. First, it is desirable that the
weight re-balancing in Step 13 due to normalization has minimal effect on the
weights of the surviving particles. We discovered through experiments that it
is good for heavy weight particles to remain heavily weighted. Therefore
$w_{new}$ should be small. Second, $w_{new}$ should be larger than
$w_{inact}$, because otherwise, the newly generated particles in a step will
be immediately deleted in the next step. Third, the purpose of setting the
value of $w_{new}$ is to give some initial weights to the new particles so
that they can participate in the weight updating in the subsequent steps. If a
new particle is fit, its weight will boost up exponentially fast; if a new
particle is unfit, it will decay exponentially fast. Therefore, the initial
weights assigned to these new particles should not significantly affect their
chance of survival and their long-term weight dynamics. Thus, as long as
$w_{new}$ is fairly small and larger than $w_{inact}$, the choice of its
actual value may not make much difference qualitatively.

More simulations are shown in Figure 16.

For the linear bandit problem, TS can also be exactly implemented by a Kalman
filter. The initial set of particles of PTS and RPTS for linear bandits are
generated uniformly at random from the unit ball in $\mathbb{R}^{K}$. That is
based on the assumption that we already know that $\theta^{*}$ is in the unit
ball before running the algorithm. In practice, such knowledge may not be
available and a common practice is to use a distribution that spreads out wide
enough so that it should cover $\theta^{*}$. For the purpose of demonstrating
the performance of PTS and RPTS here, our practice should be acceptable.

##  Appendix D Approximation of expected reward for the network slicing model

In Section 6, in step 4 of Algorithm 4, the expected reward
$\mathbb{E}_{\theta_{t}}[R(Y)|A_{t}=a,c_{t}]$ becomes
$\mathbb{E}_{\theta_{t}}[g_{c_{t,2}}(Y_{t})|a]$ for the network slicing model,
where $Y_{t}=Y_{t,1}+Y_{t,2}+Y_{t,3}$. Since $Y_{t,1},Y_{t,2},Y_{t,3}$ are
coupled through the non-linear function $g_{d}$, it is not clear if the
expectation can be exactly calculated by a closed-form expression. We propose
the following approximation. Given a random variable $Y=Y_{1}+Y_{2}+Y_{3}$,
where $Y_{i}$ is an exponentially distributed random variable with mean
$\mu_{i}$ and $Y_{i}$’s are independent. Suppose we approximate $Y$ by a
Gaussian random variable $\widetilde{Y}$ with mean
$\mu=\mu_{1}+\mu_{2}+\mu_{3}$ and variance
$\sigma^{2}=\mu_{1}^{2}+\mu_{2}^{2}+\mu_{3}^{2}$. Then

$\displaystyle\mathbb{E}[g_{d}(Y)]$
$\displaystyle\approx\mathbb{E}[g_{d}(\widetilde{Y})]$
$\displaystyle=\int_{0}^{d}\\!\frac{y}{d}\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{(y-\mu)^{2}}{2\sigma^{2}}}\,\mathrm{d}y$
$\displaystyle=\int_{-\mu}^{d-\mu}\\!\frac{1}{d}(z+\mu)\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{z^{2}}{2\sigma^{2}}}\,\mathrm{d}z\quad\left(\text{with}\;z=y-\mu\right)$
$\displaystyle=\frac{1}{d\sqrt{2\pi\sigma^{2}}}\int_{-\mu}^{d-\mu}\\!ze^{-\frac{z^{2}}{2\sigma^{2}}}\,\mathrm{d}z+\frac{\mu}{d}\int_{-\mu}^{d-\mu}\\!\frac{1}{\sqrt{2\pi\sigma^{2}}}e^{-\frac{z^{2}}{2\sigma^{2}}}\,\mathrm{d}z$
$\displaystyle=\frac{\sigma}{d\sqrt{2\pi}}\left(e^{-\frac{\mu^{2}}{2\sigma^{2}}}-e^{-\frac{(d-\mu)^{2}}{2\sigma^{2}}}\right)+\frac{\mu}{d}\left(\Phi\left(\frac{d-\mu}{\sigma}\right)-\Phi\left(-\frac{\mu}{\sigma}\right)\right)\,,$

where $\Phi(x)\triangleq\mathbb{P}(N\leq x)$ for a standard Gaussian random
variable $N$. Then

$\displaystyle\mathbb{E}_{\theta_{t}}\left[g_{c_{t,1}}(Y_{t})|a\right]\approx\frac{\sigma_{t}}{c_{t,2}\sqrt{2\pi}}\left(e^{-\frac{\mu_{t}^{2}}{2\sigma_{t}^{2}}}-e^{-\frac{(c_{t,2}-\mu_{t})^{2}}{2\sigma_{t}^{2}}}\right)+\frac{\mu_{t}}{c_{t,2}}\left(\Phi\left(\frac{c_{t,2}-\mu_{t}}{\sigma_{t}}\right)-\Phi\left(-\frac{\mu_{t}}{\sigma_{t}}\right)\right)\,,$
(23)

where $\mu_{t}=\mu_{t,1}+\mu_{t,2}+\mu_{t,3}$ and
$\sigma^{2}_{t}=\mu_{t,1}^{2}+\mu_{t,2}^{2}+\mu_{t,3}^{2}$ and
$\mu_{t,i}=c_{t,1}\theta_{t,i,a_{i},1}+\theta_{t,i,a_{i},2}$ for $i=1,2,3$.
Step 4 of Algorithm 4 can then be approximately solved by looping over all
possible $a\in[B_{1}]\times[B_{2}]\times[B_{3}]$ and find the one that
maximizes (23).