# Topology of equatorial timelike circular orbits around stationary black
holes

Shao-Wen Weia,b<EMAIL_ADDRESS>Yu-Xiao Liua,b<EMAIL_ADDRESS>aLanzhou Center for Theoretical Physics, Key Laboratory of Theoretical Physics
of Gansu Province, School of Physical Science and Technology, Lanzhou
University, Lanzhou 730000, People’s Republic of China,
bInstitute of Theoretical Physics $\&$ Research Center of Gravitation, Lanzhou
University, Lanzhou 730000, People’s Republic of China

###### Abstract

A topological approach has been successfully used to study the properties of
the light ring and the null circular orbit, in a generic black hole
background. However, for the equatorial timelike circular orbit, quite
different from the light ring case, its radius is closely dependent of the
energy and angular momentum of a test particle. This fact seems to restrict
the extension of the topological treatment to the timelike circular orbit. In
this paper, we confirm that the angular momentum does not affect the
asymptotic behavior of the constructed vector with its zero points denoting
the equatorial timelike circular orbits. As a result, a well-behaved topology
to characterize the equatorial timelike circular orbits can be constructed.
Our study shows that the total topological number of the timelike circular
orbits vanishes for a generic black hole, which is dependent of the energy of
the particle. Significantly, it reveals that if there the timelike circular
orbits exist, they always come in pairs for fixed angular momentum. Meanwhile,
the stable and unstable timelike circular orbits have positive or negative
winding number. Of particular interest is that the marginally stable circular
orbit corresponds to the bifurcation point of the zero point of the
constructed vector. Moreover, we also examine the case when the particle
energy acts as the control parameter. It is shown that there will be
topological phase transition when the value of the particle energy is one.
Below this value, the timelike circular orbits always come in pairs for fixed
energy. Otherwise, we will have one more unstable timelike circular orbit.
Such a topological phase transition actually measures whether there are
bounded orbits. We further apply the treatment to the Kerr black hole. All the
results given in a generic black hole background are exactly reproduced. These
strongly indicate that our topological approach can be generalized to the
equatorial timelike circular orbits.

Classical black hole, timelike circular orbit, topological charge

###### pacs:

04.20.-q, 04.25.-g, 04.70.Bw

##  I Introduction

Topology provides a powerful tool to understand the light rings (LRs) around
compact objects and black holes. It was shown that LRs can be well-explored by
calculating the topological charge while ignoring where the LRs are located
Cunhaa ; Cunhab . Due to the fact that the LRs are closely related to the
ringdown stage of a black hole merger and formation of black hole shadow,
further study of the topology approach may peek into the modern observations
of the LIGO/Virgo Collaboration Abbott  and EHT Collaboration Akiyama1 ;
Akiyama2 .

In 2017, Cunha, Berti, and Herdeiro Cunhaa  constructed a special vector by
using the null geodesics and found that the LRs are exactly at the zero points
of the vector. Then in the vector space, the winding number of a zero point
can provide us with some interesting topological results. Meanwhile, summing
the winding numbers of all the zero points, one can obtain the topological
charge for the LRs. Via such a topological approach, they originally proposed
that under a stationary and axisymmetric equilibrium ultracompact object, the
topological charge of the LRs always vanishes. This indicates that if the LRs
exist, they must come in pairs for such ultracompact objects. Moreover, if one
of the LRs is stable, the other one must be unstable Crispino ; Hoda .

Shortly afterwards, the study was generalized to the black hole backgrounds
Cunhab . The result suggests that for a stationary, axisymmetric,
asymptotically flat black hole spacetime in four dimensions, it possesses at
least one standard LR outside the horizon for each rotation sense. Other
related issues including the stability of angular motion were discussed in
Refs. Guo ; Junior ; Ghosh . In particular, it was found that there exists a
topological phase transition for the dilatonic Melvin spacetime and
Schwarzschilld-dilatonic-Melvin black hole at a certain value of the dilatonic
parameter JuniorYang .

Combining this with Duan’s topological current $\phi$-mapping theory Duan , we
also applied the study to the nonrotating black hole case Wei2019 . The inner
structure of the topological charge, as well as the bifurcation phenomena was
investigated. For the Dyonic black hole solution Lu , the previous study
showed that there exist several photon spheres (PSs), which are the
counterpart of LRs in spherically symmetric spacetime. This provides us with a
good opportunity to test the nontrivial topological argument. After
calculating the topological charge, we found that its value is always -1
independent of other parameters. This implies that the standard PSs are always
one more than the exotic ones. So the topological approach holds even when
multi PSs are presented. The black hole spin will also not affect the value of
the topological charge.

Inspiring by the study, we proposed a topological approach for the black hole
thermodynamics Wei2022 . Starting with the first law of black hole
thermodynamics, we constructed a vector. The critical point of the phase
transition was found to exactly locate at the zero point of the vector. By
making use of the Duan’s topological current $\phi$-mapping theory, we endowed
each critical point with a winding number. Summing the winding numbers for all
the critical points, we shall obtain the topological charge, which remains
unchanged for the system. Unexpectedly, we found that besides the conventional
critical point of negative winding number, there is another novel critical
point with positive winding number. Combining with the black hole
thermodynamics, we pointed out that only the conventional critical point acts
as an indicator of the first order black hole phase transition. Such a
topological study was also extended to the Gauss-Bonnet black hole Yerra ,
where three critical points appear. Nevertheless, the total topological charge
remains unchanged with the black hole electric charge.

Considering that the black hole solutions with the critical points are very
limited, we introduced another topological approach by the free energy of the
black hole system Weito . Each actual black hole solution corresponds to the
zero point of a constructed vector. Thus, it will be endowed with a winding
number in the vector space. Taking the inverse Hawking temperature as the
control parameter, the total topological number is universal and independent
of the black hole temperature and other parameters. By taking advantage of
such topological properties, different black hole systems are divided into
different topological classes. For example, the Schwarzschild, Reissner-
Nordström, and Reissner-Nordström-AdS black holes, respectively, have
topological number -1, 0, 1, which indicates they belong in different
topological classes. All these results indicate that the topological approach
is a promising way to study the LRs and thermodynamic properties of the black
hole systems.

Comparing with the photons, the massive particles can also go around the black
hole along the timelike circular orbits (TCOs). In certain astronomical
processes, such orbits serve as the endpoints of dynamical evolution due to
the dissipative effects Peters ; Rodrigo . In particular, there is the
innermost stable circular orbit (ISCO), below which the massive particles
shall quickly plunge into the black hole. Therefore, the stable accretion disk
can only be formed above this bound, and provide a light source to illuminate
the black hole, which may cast potential features of the black hole shadow in
the observation of EHT Akiyama1 ; Akiyama2 . For a particle falls towards the
black hole at a large radius, it will pass through a series of TCOs by losing
angular momentum, and converting its energy into radiation, and finally
reaches the ISCO before it plunges into the black hole. In the extremal Kerr
black hole background, a particle could convert its 42% energy, or achieve a
higher efficiency in compact bosonic objects Delgado . Via this process, which
leads to the possible fact that black holes could be the powerful luminosity
sources to trigger the energetic astronomical phenomena in strong gravity
regime.

On the other hand, very recently, the general properties of the TCOs have been
explored in Refs. Delgado ; Cardoso . It was found that the stable and
unstable LRs delimit the region of stable and unstable TCOs. However, the
topological approach for the TCOs remains to be established despite the great
success gained for the LRs. In the backgrounds of compact objects or black
holes, the locations of LRs are only dependent of the parameters of the
compact objects or black holes, while independent of the energy and angular
momentum of the photon. This states that the LRs are intrinsic properties of
the spacetime. However, not all the photons can surround the black hole along
these LRs. Considering that these photons come from far away from the black
hole, some of them will fall into and some of them will fly past it. Only
those with appropriate energy and angular momentum will surround the black
hole along the LRs. Moreover, if the background is fixed, the relation of the
energy and angular momentum of the photon is linear.

Now we wonder whether a massive test particle with certain energy and angular
momentum can feel the TCOs. Similar to the photon, there must be some
conditions. As is well-known, the radius of the TCO depends both on the black
hole parameters and the massive particles themselves. For a given TCO, the
relation of the energy and angular momentum of these particles orbiting along
the black hole is not linear anymore. Due to these facts, it seems that we
cannot construct a well-behaved topological approach for the TCOs as expected,
or it is impossible to obtain the corresponding topological information of the
TCOs.

However, it is worthwhile pointing out that the greatest advantage of the
topological approach is that it ignores the specific details, such as the
relation of the energy and angular momentum, as well as the locations of the
TCOs. Keeping this in mind, in order to establish a well-behaved topological
approach, we just require that, in certain parameter regions, the asymptotic
behavior of the constructed vector does not change. Motivated by this idea, we
manage to construct a topological approach for the TCOs for the first time. We
shall show that when taking the angular momentum as the control parameter, the
energy of the massive particle does not change the asymptotic behavior of the
considered vector, indicating a well-behaved topology. In such approach, the
marginally stable circular orbit (MSCO) acts as a bifurcation point of the
zero points of the vector. More importantly, we find that the TCOs with the
same angular momentum (or energy less than one) always come in pairs for the
Kerr-like black holes.

The present paper is organized as follows. In Sec. II, we start with the
geodesics and construct a vector with its zero points corresponding to the
TCOs. The asymptotic behavior of the vector is examined when we take the
angular momentum as the control parameter, which means that the TCOs evolve
with the angular momentum. In Sec. III, we calculate the global and local
topological charges for the TCOs. Meanwhile, the topological property of the
MSCO is also discussed. We then apply this approach to the Kerr black hole in
Sec. IV as an example. Adopting energy as the control parameter, we
investigate the corresponding topology in Sec. V. Finally, we summarize and
discuss our results in Sec. VI.

##  II Circular geodesics and effective potential

We consider a four-dimensional stationary, axisymmetric, asymptotically-flat
black hole spacetime possessing two Killing vectors,
$\xi^{\mu}=(\partial_{t})^{\mu}$ and $\psi^{\mu}=(\partial_{\varphi})^{\mu}$.
Moreover, we assume the metric is circular, i.e.,
$\xi^{\mu}R_{\mu}^{[\nu}\xi^{\rho}\psi^{\sigma]}=\psi^{\mu}R_{\mu}^{[\nu}\xi^{\rho}\psi^{\sigma]}=0$,
and has a north-south $\mathcal{Z}_{2}$ symmetry. As a result, the circular
orbit motion can only be on the equatorial plane. Then under these conditions,
the most generic metric reads

$\displaystyle
ds^{2}=g_{tt}(r,\theta)dt^{2}+g_{rr}(r,\theta)dr^{2}+g_{\theta\theta}(r,\theta)d\theta^{2}+g_{\varphi\varphi}d\varphi^{2}+2g_{t\varphi}dtd\varphi,$
(1)

assuming a signature (-, +, +, +). In addition to the Kerr black hole, most
known black holes can also be described by this metric. For a black hole, we
suppose that its horizon is located. at a constant (positive) radial
coordinate $r=r_{\text{h}}$. So in the exterior region of the black hole, one
should have $r_{\text{h}}<r<\infty$. We also have $g_{rr}>0$,
$g_{\theta\theta}>0$, and $g_{\varphi\varphi}>0$ outside the horizon. The
angular coordinates $\theta\in[0,\pi]$ and $\varphi\in[0,2\pi]$. According to
the $\mathcal{Z}_{2}$ symmetry, the equatorial plane is at $\theta=\pi/2$.
Although we only concern the equatorial TCOs, here we still expect to leave
$\theta$ as a free variable. Considering the condition det$(-g)>0$, one easily
obtains

$\displaystyle B(r,\theta)\equiv g_{t\varphi}^{2}-g_{tt}g_{\varphi\varphi}>0.$
(2)

At the horizon, we have $B(r_{\text{h}},\theta)=0$. Under the metric (1), the
motion of a test particle can be described by the Lagrangian

$\displaystyle\mathcal{L}=\frac{1}{2}g_{\mu\nu}\dot{x}^{\mu}\dot{x}^{\nu}=-\frac{1}{2}\mu^{2},$
(3)

where the dots denote the derivative with respect to an affine parameter, and
$\mu^{2}=$1, 0, -1 are for the timelike, null, and spacelike geodesics,
respectively. The conjugate momenta can be calculated as

$\displaystyle\pi_{\mu}=\frac{\partial\mathcal{L}}{\partial\dot{x}^{\mu}}=g_{\mu\nu}\dot{x}^{\nu}.$
(4)

Therefore, the Hamiltonian for the test particle reads

$\displaystyle\mathcal{H}$ $\displaystyle=$
$\displaystyle\pi_{\mu}\dot{x}^{\mu}-\mathcal{L}$ (5) $\displaystyle=$
$\displaystyle\frac{1}{2}\left(g_{tt}\dot{t}^{2}+g_{rr}\dot{r}^{2}+g_{\theta\theta}\dot{\theta}^{2}+g_{\varphi\varphi}\dot{\varphi}^{2}+2g_{t\varphi}\dot{t}\dot{\varphi}\right)=-\frac{1}{2}\mu^{2}.$

Rearranging it, we get

$\displaystyle
g_{rr}\dot{r}^{2}+g_{\theta\theta}\dot{\theta}^{2}+g_{tt}\dot{t}^{2}+2g_{t\varphi}\dot{t}\dot{\varphi}+g_{\varphi\varphi}\dot{\varphi}^{2}+\mu^{2}=0.$
(6)

Following the treatment for the null geodesics Grover ; Herdeiro , we denote
the kinetic term $\mathcal{K}$ and potential term $\mathcal{V}$ as

$\displaystyle\mathcal{K}$ $\displaystyle=$ $\displaystyle
g_{rr}\dot{r}^{2}+g_{\theta\theta}\dot{\theta}^{2},$ (7)
$\displaystyle\mathcal{V}$ $\displaystyle=$ $\displaystyle
g_{tt}\dot{t}^{2}+2g_{t\varphi}\dot{t}\dot{\varphi}+g_{\varphi\varphi}\dot{\varphi}^{2}+\mu^{2}.$
(8)

Therefore, the motion of the test particle (6) turns to

$\displaystyle\mathcal{K}+\mathcal{V}=0.$ (9)

Note that the kinetic term $\mathcal{K}\geq 0$, and the inequality is only
saturated at $\dot{r}=\dot{\theta}=0$. Then the motion of the particle is
completely governed by the effective potential $\mathcal{V}$. On the other
hand, these two Killing vectors, $\xi^{\mu}$ and $\psi^{\mu}$, are related to
two constants, the energy and orbital angular momentum of the test particle,
along each geodesic

$\displaystyle-E$ $\displaystyle=$ $\displaystyle
g_{\mu\nu}u^{\mu}\xi^{\nu}=g_{tt}\dot{t}+g_{t\varphi}\dot{\varphi},$ (10)
$\displaystyle l$ $\displaystyle=$ $\displaystyle
g_{\mu\nu}u^{\mu}\psi^{\nu}=g_{t\varphi}\dot{t}+g_{\varphi\varphi}\dot{\varphi},$
(11)

where $u^{\mu}$ denotes the tangent vector of the geodesics. Solving the above
two equations, one can obtain

$\displaystyle\dot{t}$ $\displaystyle=$
$\displaystyle\frac{1}{B}\left(Eg_{\varphi\varphi}+lg_{t\varphi}\right),$ (12)
$\displaystyle\dot{\varphi}$ $\displaystyle=$
$\displaystyle-\frac{1}{B}\left(Eg_{t\varphi}+lg_{tt}\right).$ (13)

Plugging them into the potential (8), one has

$\displaystyle\mathcal{V}=-\frac{1}{B}\left(E^{2}g_{\varphi\varphi}+2Elg_{t\varphi}+l^{2}g_{tt}\right)+\mu^{2}.$
(14)

When $\mu^{2}=0$, this exactly reduces to the potential of a photon Grover ;
Herdeiro . Similarly, we can reexpress this potential as

$\displaystyle\mathcal{V}=-\frac{l^{2}g_{\varphi\varphi}}{B}\left(\frac{E}{l}-H_{+}\right)\left(\frac{E}{l}-H_{-}\right)+\mu^{2},$
(15)

where $H_{\pm}$ are only dependent of the black hole metric and are given by

$\displaystyle H_{\pm}=\frac{-g_{t\varphi}\pm\sqrt{B}}{g_{\varphi\varphi}}.$
(16)

For a photon with $\mu^{2}=0$, it is quite clear that the LR related with
$\mathcal{V}=\partial_{r}\mathcal{V}=0$ satisfies Grover ; Herdeiro

$\displaystyle\frac{E}{l}$ $\displaystyle=$ $\displaystyle H_{\pm},$ (17)
$\displaystyle\frac{\partial H_{\pm}}{\partial r}$ $\displaystyle=$
$\displaystyle 0.$ (18)

This is an amazing result that the radius $r_{LR}$ of the LR (or the photon
sphere for the static spherically black hole) of the black hole can be solved
by Eq. (18) and is independent the properties of the photon, indicating that
it is an intrinsic structure of the black hole just like the horizon. Then
plugging it into (17), one shall get $E/l=H_{\pm}(r_{LR})$ for these photons
orbiting along the LR. In Ref. Cunhab , by defining a new vector
$v=((\partial_{r}H_{\pm})/\sqrt{g_{rr}},(\partial_{\theta}H_{\pm})/\sqrt{g_{\theta\theta}})$,
one finds the location of the LR is exactly at the zero point of $v$ from
(18). As a result, by calculating the winding number in the vector space, we
can obtain some properties of the LR. In particular, considering all the
regions of ($r$, $\theta$), the global property of the LR in a black hole
background will be obtained by examining the asymptotic behaviors of $v$ near
the boundary. This treatment mainly takes advantage of the result that
$H_{\pm}$ are independent of the energy and angular momentum of the photon.
However, due to the presence of the extra term $\mu^{2}$ in (15), the case
becomes unclear.

In order to establish the topology for the TCO, we reexpress the effective
potential (14) in the following form

$\displaystyle\mathcal{V}=-\frac{g_{\varphi\varphi}}{B}(E-e_{1})(E-e_{2}),$
(19)

where

$\displaystyle
e_{1}=\frac{-lg_{t\varphi}+\sqrt{B}\sqrt{l^{2}+\mu^{2}g_{\varphi\varphi}}}{g\varphi\varphi},$
(20) $\displaystyle
e_{2}=\frac{-lg_{t\varphi}-\sqrt{B}\sqrt{l^{2}+\mu^{2}g_{\varphi\varphi}}}{g\varphi\varphi}.$
(21)

Obviously, $\mathcal{V}=0$ leads to the energy of the test particle
$E=e_{1,2}$. Considering the condition (2) and $g_{\varphi\varphi}>0$, it is
easy to obtain $e_{2}<g_{t\varphi}(|l|-l)/g_{\varphi\varphi}$. For $l\geq 0$,
it is easy to obtain $e_{2}<0$, and for $l<0$,
$e_{2}<-2lg_{t\varphi}/g_{\varphi\varphi}<0$ due to negative $g_{t\varphi}$.
Therefore, it is notable that for any value of $l$, $e_{2}$ takes a negative
value. As a result, we abandon it. Then, the conditions
$\mathcal{V}=\partial_{r}\mathcal{V}=0$ for determining the TCO become

$\displaystyle E$ $\displaystyle=$ $\displaystyle e_{1},$ (22)
$\displaystyle\partial_{r}e_{1}$ $\displaystyle=$ $\displaystyle 0.$ (23)

Quite different from the LR, $\partial_{r}e_{1}$ here not only depends on $r$,
but also on the angular momentum $l$ of the test particle, so it seems that we
cannot construct the topology for the TCO as expected.

However, we would like to propose a new viewpoint while is still similar to
that for the photon. For each given angular momentum $l$, we can obtain the
radius of the TCO by solving (23). Whether the equation has solution or not,
the angular momentum $l$ of the test particle does not affect our treatment.
After obtaining the radius of the TCO, we can substitute it into (22), and
then the energy of the test particle will be obtained.

In order to give a global topology, we only require that the values of the
angular momentum do not change the asymptotic behavior of $\partial_{r}e_{1}$
at the boundary of the ($r$, $\theta$) plane. Aiming at it, we define a new
vector $\phi$=($\phi^{r}$, $\phi^{\theta}$) by analog to Ref. Cunhab  for the
photon,

$\displaystyle\phi^{r}=\frac{\partial_{r}e_{1}}{\sqrt{g_{rr}}},\quad\phi^{\theta}=\frac{\partial_{\theta}e_{1}}{\sqrt{g_{\theta\theta}}},$
(24)

in a flat vector space. In order to be consistent with the notations of Ref.
Duan , we write these indices with up ones. Obviously, $\phi=0$ corresponds to
the TCO and $\theta=\pi/2$ as expected. Next, we turn to the asymptotic
behavior of $\phi$ at the boundary of the ($r$, $\theta$) plane; this shall
further confirm that the construction of the topology for the TCO is
meaningful.

For our case, we should examine four asymptotic behaviors of $\phi$ near the
horizon $r\rightarrow r_{\text{h}}$, infinity $r\rightarrow\infty$, and
$\theta\rightarrow 0$, $\pi$.

i) _Axis limits_. Here we shall consider the axis limits of the vector $\phi$.
In order to clearly show them, we here follow Ref. Cunhab  and give a detailed
calculation. Defining a local coordinate $\rho=\sqrt{g_{\varphi\varphi}}$, one
can find when $\theta\rightarrow 0$ or $\pi$, $d\rho/d\theta$ is positive or
negative. Since the metric functions are independent of the particles, one has
the same result as that of Ref. Carot  by expanding the metric function at
small $\rho$

$\displaystyle g_{tt}\sim g_{tt}^{0}+\mathcal{O}(\rho),\quad g_{\rho\rho}\sim
g_{\rho\rho}^{0}+\mathcal{O}(\rho),$ (25) $\displaystyle
g_{\varphi\varphi}=\rho^{2},\quad g_{t\varphi}\sim
b_{0}\rho^{n}+\mathcal{O}(\rho^{n+1}),$ (26)

where $g_{tt}^{0}$, $g_{\rho\rho}^{0}$, $b_{0}$ are constants, and $n\in N$.
Assuming $C^{2}$ smoothness and regularity close to the axis, one shall have
$n\geq 2$ Carot . Plugging them into Eq. (20), we obtain

$\displaystyle
e_{1}\sim\frac{\sqrt{-g_{tt}^{0}l^{2}\rho^{2}+\sqrt{-g_{tt}^{0}}\mu^{2}\rho^{3}}}{\rho^{2}}.$
(27)

At zero order in $\rho$, one has $g_{\rho\rho}d\rho^{2}\sim
g_{\theta\theta}d\theta^{2}$, and thus

$\displaystyle\phi^{\theta}\sim\text{sign}\left(\frac{d\rho}{d\theta}\right)\partial_{\rho}e_{1}\sim-\text{sign}\left(\frac{d\rho}{d\theta}\right)\frac{|l|\sqrt{-g_{tt}^{0}}}{\rho^{2}},$
$\displaystyle\phi^{\theta}\sim\text{sign}\left(\frac{d\rho}{d\theta}\right)\partial_{\rho}e_{1}\sim-\text{sign}\left(\frac{d\rho}{d\theta}\right)\frac{(-g_{tt}^{0})^{\frac{1}{4}}\mu^{2}}{2\rho^{\frac{3}{2}}},$
(28)

for nonvanishing and vanishing $l$, respectively. Therefore, it is clear that
$\phi^{\theta}$ diverges at $\rho\to 0$ and the angular momentum do not affect
the behavior of $\phi^{\theta}$ at small $\rho$. Further combining with
$\phi^{r}\sim\rho^{-1}$, we have $\phi^{\theta}\gg\phi^{r}$. This indicates
that the direction of the vector is vertical near $\theta$=0 and $\pi$.
Moreover, considering the fact that $d\rho/d\theta$ is positive and negative
as $\theta\rightarrow$0 and $\pi$, we find the vector is outward at $\theta$=0
and $\pi$.

ii) _Horizon limit_. Here we would like to consider the behavior of the vector
$\phi$ when $r\to r_{h}$. As discussed in Ref. Medved , we can set a local
radial coordinate $x$ such that $g_{xx}=1$, and $x=0$ at the black hole
horizon. Near the horizon, one can obtain these metric functions Medved

$\displaystyle\omega\simeq\omega_{H}+\mathcal{O}(x^{2}),\quad
g_{\varphi\varphi}\simeq g_{\varphi\varphi}^{H}+\mathcal{O}(x^{2}),$ (29)
$\displaystyle g_{t\varphi}\simeq-\omega_{H}g_{\varphi\varphi}^{H},\quad
B\simeq g_{\varphi\varphi}^{H}\kappa^{2}x^{2},\quad
g_{tt}\simeq\omega_{H}^{2}g_{\varphi\varphi}^{H},$ (30)

where $\omega=-g_{t\varphi}/g_{\varphi\varphi}$ and $\kappa$ denotes the
surface gravity of the black hole. $g_{\varphi\varphi}^{H}$ and $\omega_{H}$
are zero-order terms of the expansions in $x$. Then

$\displaystyle\partial_{x}e_{1}\simeq\kappa\sqrt{\frac{l^{2}}{g_{\varphi\varphi}^{H}}+\mu^{2}}+\mathcal{O}(x).$
(31)

By using $(1/\sqrt{g_{xx}})(\partial/\partial
x)=(1/\sqrt{g_{rr}})(\partial/\partial r)$, one shall obtain

$\displaystyle\phi^{r}=\frac{\partial_{r}e_{1}}{\sqrt{g_{rr}}}\simeq\kappa\sqrt{\frac{l^{2}}{g_{\varphi\varphi}^{H}}+\mu^{2}}.$
(32)

For a nonextremal black hole with $\kappa>0$, we see that $\phi^{r}$ takes a
positive value. In order to determine the direction of the vector $\phi$ at
the horizon, one needs to calculate $\phi^{\theta}$. However, no matter
whether $\phi^{\theta}$ is positive or negative, we can claim that the
direction of the vector $\phi$ must be to the right, with at most an
inclination.

iii) _Asymptotic limit_. Here let us consider the asymptotic behavior of the
vector at infinity. Here we require the spacetime tends to a flat one in
standard spherical coordinates, namely

$\displaystyle g_{t\varphi}$ $\displaystyle\simeq$
$\displaystyle-\mathcal{O}\left(\frac{1}{r}\right),\quad
g_{\varphi\varphi}\simeq r^{2}\sin^{2}\theta+\mathcal{O}\left(1\right),$
$\displaystyle g_{tt}$ $\displaystyle\simeq$
$\displaystyle-1+\frac{2M}{r}+\mathcal{O}\left(\frac{1}{r^{2}}\right),\quad
g_{rr}\simeq 1+\frac{2M}{r}++\mathcal{O}\left(\frac{1}{r^{2}}\right).$ (33)

Parameter $M$ is a positive constant and can be treated as the black hole mass
in the asymptotically flat spacetime. Therefore, we have

$\displaystyle\phi^{r}\simeq\frac{M\sqrt{\mu^{2}}}{r^{2}}+\mathcal{O}\left(\frac{1}{r^{3}}\right),$
(34)

on the equatorial plane by taking $\theta=\pi/2$, which is small but a
positive value. As a result, by ignoring the specific value of
$\phi^{\theta}$, the direction of the vector $\phi$ is to the right similar to
that in the horizon limit.

In summary, we see that the angular momentum $l$ of the test massive particle
do not affect the behavior of the vector $\phi$ at the boundary of the ($r$,
$\theta$) plane in our above approach. This result strongly implies that we
can construct the topology for the TCOs of the massive particle in a black
hole background.

##  III Global and local topologies

Based on the asymptotic behavior of $\phi$ at the boundary of the ($r$,
$\theta$) plane, we here investigate the global and local topological
properties for the TCOs. Moreover, we also examine the topology for the MSCO
and ISCO.

###  III.1 Global property

As we have shown that the TCOs locate exactly at the zero points of the vector
$\phi$. Hence, we can calculate the topological number $W$ defined as Duan

$\displaystyle W=\int_{\Sigma}j^{0}d^{2}x$ (35)

for a giving region $\Sigma$ in the ($r$, $\theta$) plane for the TCOs. Here
$j^{0}$ is the zeroth component of the topological current

$\displaystyle
j^{\mu}=\frac{1}{2\pi}\epsilon^{\mu\nu\rho}\epsilon_{ab}\frac{\partial
n^{a}}{\partial x^{\nu}}\frac{\partial n^{b}}{x^{\rho}},$ (36)

where $x^{\mu}=(t,r,\theta)$, and the unit vector $n^{a}$=($n^{r}$,
$n^{\theta}$)=($\phi^{r}/|\phi|$, $\phi^{\theta}/|\phi|$). Here, $x^{0}=t$
denotes the time control parameter. One can study the evolution of the zero
points with that parameter as expected. On the other hand, other parameters of
the black hole or particle can also be chosen as the control parameter. For a
good definition, we require that the zero points cannot reach the boundary of
the ($r$, $\theta$) plane at finite value of the control parameter. In the
following study, we mainly treat the angular momentum $l$ of the particle as a
control parameter. The brief discussion on choosing energy will also be given.

After some calculations, one arrives at Duan

$\displaystyle W=\sum_{n=1}^{N}w_{n},$ (37)

which is the sum of the winding number $w$ for all the zero points on the
region $\Sigma$. On the other hand, we can also obtain the topological charge
for a region $\Sigma$ by counting the number of the loops that the vector
makes in the vector $\phi$ space when $x^{\mu}$ moves along the closed curve
$\partial\Sigma$ in a counterclockwise direction,

$\displaystyle
W=\frac{1}{2\pi}\oint_{C}d\Omega=\frac{1}{2\pi}\oint_{C}\epsilon_{ab}n^{a}dn^{b},$
(38)

where $\Omega$ denotes the change of the direction of the vector.

Figure 1: Representation of the contour $C=\partial\Sigma=L_{1}\cup L_{2}$
(which encloses all the parameter region $\Sigma$) on the ($r$, $\theta$)
plane. $L_{1}$ and $L_{2}$ are the left and right parts of the contour marked
in blue and red colors, respectively. Black arrows indicate the approximate
directions of the vector $\phi$ in the boundaries.

Considering that the curve $C=\partial\Sigma$ encloses the entire parameter
region, one will obtain the global property of the topology. Here we can
divide the closed curve into two parts $L_{1}$ and $L_{2}$ such that
$C=L_{1}\cup L_{2}$. As shown in Fig. 1, we can see that the direction of the
vector are outward at $\theta=0$ and $\pi$, and is to the right at
$r=r_{\text{h}}$ and $\infty$ by ignoring the possible inclination. According
to the asymptotic behaviors of $\phi$, we easily get

$\displaystyle W$ $\displaystyle=$
$\displaystyle\frac{1}{2\pi}\int_{L_{1}}d\Omega+\frac{1}{2\pi}\int_{L_{2}}d\Omega$
(39) $\displaystyle=$
$\displaystyle\frac{1}{2\pi}\times(-\pi)+\frac{1}{2\pi}\times(\pi)$
$\displaystyle=$ $\displaystyle 0.$

This states that the total topological number vanishes for the TCOs.
Significantly, it is universal result and is independent of the angular
momentum of the massive test particle, as well as the black hole spin. This is
an amazing result, and based on which, we can conclude that the TCOs with the
fixed angular momentum always come in pairs for the four-dimensional
stationary, axisymmetric, asymptotically flat black hole. This property is
quite different from that of the LR, whose topological number is -1 indicating
that at least one one LR exists.

Figure 2: Different patterns for the direction of the vector $\phi$ at
$\theta$=0 and $\pi$.

As shown above, the total topological number for the TCOs vanishes where the
direction of the vector $\phi$ is outwards at $\theta$=0 and $\pi$. However,
even under the $\mathcal{Z}_{2}$ symmetry, one may wonder whether the
topological number still vanishes if the direction of the vector $\phi$ is
inwards at $\theta$=0 and $\pi$, see Fig. 2. After considering the change of
the direction of the vector $\phi$, we still obtain $W=0$. On the other hand,
it is valuable to examine the total topological number if we abandon the
$\mathcal{Z}_{2}$ symmetry. We further exhibit these two possible cases in
Figs. 2 and 2, where the $\mathcal{Z}_{2}$ symmetry is broken as expected.
Surprisingly, we still have

$\displaystyle W=0,$ (40)

for both the cases without $\mathcal{Z}_{2}$ symmetry.

In short, we find that the result $W=0$ always holds for the case with or
without $\mathcal{Z}_{2}$ symmetry.

###  III.2 Local property and stability

As shown above, we confirm that the global topological number of the TCOs is
zero. So they either disappear or appear in pairs. However, quite different
from the LR which is the intrinsic structure of the spacetime, TCOs are
closely dependent of the test particle. For some certain values of the energy
and angular momentum, the TCOs may be present, while for other values, they
may be absent. In this subsection, we focus on the local topological charge
and the stability for the TCOs.

For each TCO at a given radius $r_{\text{t}}$, the energy and angular momentum
are not independent of each other, and they must satisfy the conditions (22)
and (23). From (23), one can solve the angular momentum as Delgado

$\displaystyle
l_{\pm}=\frac{g_{t\varphi}+g_{\varphi\varphi}\Omega_{\pm}}{\sqrt{\beta_{\pm}}}\bigg{|}_{r_{\text{t}}},$
(41)

where

$\displaystyle\Omega_{\pm}$ $\displaystyle=$
$\displaystyle\frac{-g_{t\varphi}^{\prime}\pm\sqrt{C}}{g_{\varphi\varphi}^{\prime}},$
(42) $\displaystyle\beta_{\pm}$ $\displaystyle=$ $\displaystyle-
g_{tt}-2g_{t\varphi}\Omega_{\pm}-g_{\varphi\varphi}\Omega_{\pm}^{2},$ (43)
$\displaystyle C$ $\displaystyle=$
$\displaystyle(g_{t\varphi}^{\prime})^{2}-g_{tt}^{\prime}g_{\varphi\varphi}^{\prime},$
(44)

with the prime denoting the derivative with respect to the radial coordinate.
Plugging the angular momentum $l$ into (22), one obtains the energy

$\displaystyle
E_{\pm}=-\frac{g_{tt}+g_{t\varphi}\Omega_{\pm}}{\sqrt{\beta_{\pm}}}\bigg{|}_{r_{\text{t}}}.$
(45)

By making use of (19), we have

$\displaystyle\frac{\partial^{2}\mathcal{V}}{\partial
r^{2}}=\frac{(e_{1}-e_{2})g_{\varphi\varphi}\sqrt{g_{rr}}}{B}\frac{\partial\phi^{r}}{\partial
r},$ (46)

for the TCOs. Since $e_{2}<0$, the coefficient is positive. Therefore, it is
easy to know

$\displaystyle\frac{\partial^{2}\mathcal{V}}{\partial
r^{2}}\bigg{|}_{r_{\text{t}}}>0(<0)\Rightarrow\frac{\partial\phi^{r}}{\partial
r}\bigg{|}_{r_{\text{t}}}>0(<0).$ (47)

Note that the stable and unstable TCOs correspond to positive and negative
$\partial^{2}\mathcal{V}/\partial r^{2}$, respectively. Via (47), we can also
examine the radial stability of the TCOs by using $\phi^{r}$. In Fig. 3 we
show an example that there are two TCOs located at $r_{1}$ and $r_{2}$
(suppose $r_{1}<r_{2}$), respectively. Via the asymptotic behaviors of
$\phi^{r}$ at horizon and infinity, we know $\phi^{r}$ is positive in the
regions $(r_{\text{h}},r_{1})\cup(r_{2},+\infty)$, and negative in
$(r_{1},r_{2})$. So near $r=r_{1}$, $\partial\phi^{r}/\partial r$ is negative,
while near $r=r_{2}$, $\partial\phi^{r}/\partial r$ is positive. This
indicates that the TCO is unstable at $r_{1}$, while stable at $r_{2}$.

Here, we would like to briefly discuss the local stability in the angular
direction. Here we only consider the equatorial TCOs with $\theta=\pi/2$. As
we have shown in Sec. II, the presence of $l$ does not change the axis limit
of the vector at $\theta$=0 and $\pi$. Moreover, we can obtain the result that
the local stability in the vertical direction of the TCOs is the same as that
of the LRs. Then according to Refs. Cunhab ; Guo , it is easy to find that the
TCOs are stable in the vertical direction.

Figure 3: Representation of the direction of the vector $\phi$ along the
contours enclosing two zero points at $r_{1}$ and $r_{2}$.

After addressing the stability of the TCOs, we expect to examine the value of
the winding number for these two TCOs. Considering that the component
$\phi^{r}$ of vector $\phi$ or $n^{r}$ of the unit vector $n$ changes its sign
at $r_{1}$ and $r_{2}$, we also mark the direction of the vector in Fig. 3.
Counting the change of the direction of the vector, one can easily get the
winding number $w$=-1 and 1 for the zero points located at $r_{1}$ and
$r_{2}$, respectively. So we can conclude that the unstable and stable TCOs
correspond to negative and positive winding numbers. Nevertheless, the sum of
the winding numbers vanishes.

###  III.3 MSCO and bifurcation point

We have shown above that positive or negative winding numbers correspond to
stable or unstable TCO with $\partial^{2}\mathcal{V}/\partial r^{2}>0$ or
$<0$. Between these two cases, there is a critical case with
$\partial^{2}\mathcal{V}/\partial r^{2}=0$, which exactly corresponds to the
MSCO or ISCO for the Kerr-like black holes.

As noted in Refs. Delgado ; Raduba ; Raduad , MSCO and ISCO are two different
concepts for (hairy) black hole or generic ultracompact objects. The MSCO
denotes the stable TCO with the smallest radius, which can be continuously
connected to spatial infinity by a set of stable TCOs. While ISCO denotes the
stable TCO with the smallest radius among all the TCOs, which can be
determined by $C=0$ (44). For the Schwarzschild and Kerr black holes, the MSCO
meets the ISCO; however, for other rotating black holes, such as the boson
stars, Proca stars, and the hairy black holes, they are different. In order to
avoid confusing, we here adopt the treatment of Ref. Delgado .

In general, the MSCO is determined by

$\displaystyle\mathcal{V}=0,\quad\frac{\partial\mathcal{V}}{\partial
r}=0,\quad\frac{\partial^{2}\mathcal{V}}{\partial r^{2}}=0.$ (48)

The first condition is naturally satisfied with $E=e_{1}$, which is dependent
of the angular momentum $l$ and radius $r_{\text{t}}$ of the TCOs. Meanwhile
the radius $r_{\text{t}}$ can be obtained by solving the second equation. As
expected, we denote its solution as $r=r_{\text{t}}(l)$. For different values
of $l$, the TCOs have different orbit radii, so we can treat the angular
momentum as a control parameter. After a simple calculation, the third
condition $\partial^{2}\mathcal{V}/\partial r^{2}=0$ or
$\partial\phi^{r}/\partial r=0$ turns to

$\displaystyle\frac{dl}{dr_{\text{t}}}\bigg{|}_{(l^{*},r_{\text{t}}^{*})}=0.$
(49)

According to Duan’s topological current theory, the point
$(l^{*},r_{\text{t}}^{*})$ related by the MSCO is a bifurcation point. Near
the point, we can expand the angular momentum as

$\displaystyle
l-l^{*}=\frac{1}{2}\frac{d^{2}l}{dr_{\text{t}}^{2}}\bigg{|}_{(l^{*},r_{\text{t}}^{*})}(r_{\text{t}}-r_{\text{t}}^{*})^{2}+\mathcal{O}(r_{\text{t}}-r_{\text{t}}^{*})^{3}.$
(50)

Via this equation, one can obtain two branch solutions in the $r_{\text{t}}-l$
plane. When
$\frac{d^{2}l}{dr_{\text{t}}^{2}}\big{|}_{(l^{*},r_{\text{t}}^{*})}>0$, these
two branch solutions are present for $l>l^{*}$, and otherwise, they are
present for $l<l^{*}$. According to this property, the bifurcation point can
be, respectively, divided into the generated and annihilated points Fu .

Here, let us focus on the topological number of the MSCO. To clearly show it,
we take the BP2 branch described by the blue curve in Fig. 7 as an example.
For the large angular momentum, there are two TCOs, while no TCO for small
one. This implies that BP2 is a generated point. As we have shown above, for
the fixed angular momentum, the TCOs with large and small radius are stable
and unstable, and which have positive winding number $w$=1 and -1,
respectively. Therefore, the upper branch has positive winding number, while
the lower branch has negative winding number. Obviously, we can see that the
branches with positive slope and negative slope, or equivalently with
$\frac{\partial^{2}\mathcal{V}}{\partial r^{2}}>0$ and $<0$ have $w$=1 and -1,
respectively. Starting with a large angular momentum, the upper and lower
branches have $\frac{\partial^{2}\mathcal{V}}{\partial r^{2}}>0$ and $<0$.
This property holds with the decrease of the angular momentum. However, when
the point BP2 is reached, such pattern changes. These two branches merge and
only one MSCO with $\frac{\partial^{2}\mathcal{V}}{\partial r^{2}}=0$ leaves.
So one can regard the MSCO as one degenerated TCO. For any angular momentum
slightly larger than that of BP2, the MSCO will be split into two TCOs, one of
which have $w$=1, while other one has $w$=-1. Considering that
$\frac{\partial^{2}\mathcal{V}}{\partial r^{2}}>0$ and $<0$, respectively,
correspond to $w$=+1 and -1. The MSCO with
$\frac{\partial^{2}\mathcal{V}}{\partial r^{2}}=0$ will not take $w$=+1 or -1.
Therefore, The MSCO must have $w$=0. This result is exactly consistent with
the result of Ref. Duan . Thus, we can say that the topological charge of the
MSCO vanishes. It can also be further confirmed in Fig. 6 and Fig. 6 by
calculating the change of the vector direction along the constructed closed
loop.

##  IV Timelike circular orbits and topology in Kerr black hole

In the previous section, we constructed the topology for the TCOs in a generic
black hole background. In this section, we manage to apply it to the Kerr
black hole and to confirm whether the above study is applicable. This shall
uncover the detailed properties for a specific example.

###  IV.1 Effective potential

The metric of the Kerr black hole is given by

$\displaystyle ds^{2}=-\frac{\Delta}{\rho^{2}}\bigg{(}dt-a\sin^{2}\theta
d\varphi\bigg{)}^{2}+\frac{\rho^{2}}{\Delta}dr^{2}+\rho^{2}d\theta^{2}$
$\displaystyle+\frac{\sin^{2}\theta}{\rho^{2}}\bigg{(}adt-(r^{2}+a^{2})d\varphi\bigg{)}^{2},$
(51)

where the metric functions read

$\displaystyle\rho^{2}$ $\displaystyle=$ $\displaystyle
r^{2}+a^{2}\cos^{2}\theta,$ (52) $\displaystyle\Delta$ $\displaystyle=$
$\displaystyle r^{2}-2Mr+a^{2}.$ (53)

Solving $\Delta=0$, one can obtain the radii of the black hole horizons

$\displaystyle r_{\pm}=M\pm\sqrt{M^{2}-a^{2}}.$ (54)

Obviously, for a black hole, its spin must be $a/M\in(0,1)$. From (19), the
effective potential $\mathcal{V}$ can be obtained, and $e_{1,2}$ are given by

$\displaystyle
e_{1,2}=\frac{2alMr\pm\rho\csc\theta\sqrt{\Delta\left(\left(a^{2}+r^{2}\right)^{2}\sin^{2}\theta+l^{2}\rho^{2}-a^{2}\Delta\sin^{4}\theta\right)}}{\left(a^{2}+r^{2}\right)^{2}-a^{2}\Delta\sin^{2}\theta},$
(55)

where we have taken $\mu^{2}=1$. In the equatorial plane, they reduce to

$\displaystyle e_{1,2}=\frac{2alMr\pm
r\sqrt{\Delta\left(\left(a^{2}+r^{2}\right)^{2}+l^{2}r^{2}-a^{2}\Delta\right)}}{\left(a^{2}+r^{2}\right)^{2}-a^{2}\Delta}.$
(56)

Figure 4: Behavior of the effective potential $\mathcal{V}$ with $a$=0.98 and
$M=1$. (a) $l$=-4.5. $E$=0.9650, 0.9691, 0.9750, 0.9776, and 0.9850 from top
to bottom. (b) $l$=-4.2209. $E$=0.9450, 0.9500, 0.9620, 0.9700, and 0.9800
from top to bottom. (c) $l$=-4. $E$=0.9450, 0.9500, 0.9620, 0.9700, and 0.9800
from top to bottom. (d) $l$=1.5. $E$=0.65, 0.70, 0.75, 0.80, and 0.85 from top
to bottom. (e) $l$=1.6827. $E$=0.75, 0.76, 0.7661, 0.77, and 0.78 from top to
bottom. (f) $l$=1.9. $E$=0.81, 0.8210, 0.83, 0.8474, and 0.86 from top to
bottom. Black dots denote these points satisfying
$\mathcal{V}=\partial_{r}\mathcal{V}$=0.

In order to show the characteristic behavior of the potential, we plot it in
Fig. 4 for the black hole with $a$=0.98 and $M=1$. When $l$=-4.5, we observe
from Fig. 4 that there are two extremal points for fixed energy $E$. When
$E$=0.9691 and 0.9776, there exists one point with
$\mathcal{V}=\partial_{r}\mathcal{V}$=0 for each curve. It is obvious that
these points are just the locations of TCOs, where the surrounding particles
have zero radial velocity. So there are two TCOs for $l$=-4.5. Slightly
increasing the angular momentum such that $l$=-4.2209, we can find from Fig. 4
that there are two extremal points for the small value of the energy. With the
increase of the energy, these two extremal points coincide at the black dot
with energy $E$=0.9620, and an extra condition $\partial_{r,r}\mathcal{V}=0$
is satisfied. Actually, the black dot denotes the MSCO as expected. In
particular, it is worthwhile pointing out that such MSCO exactly meets the
ISCO for the Kerr black hole. When the energy is beyond 0.9620, no extremal
point will be present. Taking $l$=-4, and 1.5, we describe the effective
potentials in Figs. 4 and 4, respectively. Although the angular momenta have
different signs, their patterns are quite similar. Even varying the energy, we
see that there is no extremal point, which strongly indicates that TCOs do not
exist for these cases. When $l$=1.6827 shown in Fig. 4, the effective
potential shows a similar behavior as that of Fig. 4. The MSCO is present for
the certain values of $a$, $l$, and $E$. For $l$=1.7, two TCOs marked with the
black dots are shown in Fig. 4. These disclose the characteristic behaviors of
the effective potential. In summary, by varying the angular momentum, we find
that there can be two TCOs, one MSCO, or no TCO.

###  IV.2 Vector and asymptotic behaviors

As shown above, the TCOs are related with the zero points of the vector. Here
we would like to examine the asymptotic behavior of the vector in the
background of the Kerr black hole and to see whether it is consistence with
the generic case shown in Sec. II.

The vector can be obtained via (24), which is in a complicated form and we
will not show them. Expanding the vector near $\theta=0$, we have

$\displaystyle\phi^{r}(\theta\rightarrow 0)$ $\displaystyle=$
$\displaystyle-\frac{l\left(a^{2}(M+r)+r^{2}(r-3M)\right)}{\left(a^{2}+r^{2}\right)^{\frac{5}{2}}}\theta^{-1}+\frac{2a\sqrt{{\Delta}}lM\left(a^{2}-3r^{2}\right)}{\left(a^{2}+r^{2}\right)^{\frac{7}{2}}}+\mathcal{O}(\theta),$
(57) $\displaystyle\phi^{\theta}(\theta\rightarrow 0)$ $\displaystyle=$
$\displaystyle-\frac{l\sqrt{\Delta}}{\left(a^{2}+r^{2}\right)^{\frac{3}{2}}}\theta^{-2}$
(58)
$\displaystyle+\frac{\sqrt{\Delta}\left(3a^{6}+a^{4}\left(9r^{2}-2l^{2}\right)+a^{2}r\left(9r^{3}-l^{2}(12M+r)\right)+r^{4}\left(l^{2}+3r^{2}\right)\right)}{6l\left(a^{2}+r^{2}\right)^{\frac{7}{2}}}+\mathcal{O}(\theta),$

and near $\theta=\pi$, we have

$\displaystyle\phi^{r}(\theta\rightarrow\pi)$ $\displaystyle=$
$\displaystyle\frac{l\left(a^{2}(M+r)+r^{2}(r-3M)\right)}{\left(a^{2}+r^{2}\right)^{\frac{5}{2}}}(\theta-\pi)^{-1}+\frac{2a\sqrt{{\Delta}}lM\left(a^{2}-3r^{2}\right)}{\left(a^{2}+r^{2}\right)^{\frac{7}{2}}}+\mathcal{O}(\theta-\pi),$
(59) $\displaystyle\phi^{\theta}(\theta\rightarrow\pi)$ $\displaystyle=$
$\displaystyle\frac{l\sqrt{\Delta}}{\left(a^{2}+r^{2}\right)^{\frac{3}{2}}}(\theta-\pi)^{-2}$
(60)
$\displaystyle-\frac{\sqrt{\Delta}\left(3a^{6}+a^{4}\left(9r^{2}-2l^{2}\right)+a^{2}r\left(9r^{3}-l^{2}(12M+r)\right)+r^{4}\left(l^{2}+3r^{2}\right)\right)}{6l\left(a^{2}+r^{2}\right)^{\frac{7}{2}}}+\mathcal{O}(\theta-\pi).$

Therefore, at $\theta=0$ and $\theta=\pi$, the direction
$\arg\phi=\arctan(\phi^{\theta}/\phi^{r})$ of the vector is

$\displaystyle\arg\phi(\theta\rightarrow 0)$ $\displaystyle\sim$
$\displaystyle\arctan\left(-\frac{l\sqrt{\Delta}}{\left(a^{2}+r^{2}\right)^{\frac{3}{2}}}\theta^{-2}\right)\sim\arctan\left(-\theta^{-2}\right)=-\frac{\pi}{2},$
(61) $\displaystyle\arg\phi(\theta\rightarrow\pi)$ $\displaystyle\sim$
$\displaystyle\arctan\left(\frac{l\sqrt{\Delta}}{\left(a^{2}+r^{2}\right)^{\frac{3}{2}}}\theta^{-2}\right)\sim\arctan\left(\theta^{-2}\right)=\frac{\pi}{2}.$
(62)

As a result, the direction of the vector is up at $\theta=\pi$ and is down at
$\theta=0$. On the equatorial plane, we can also expand the vector at large
$r$, which gives

$\displaystyle\phi^{r}(\theta=\frac{\pi}{2},r\rightarrow\infty)=\frac{M}{r^{2}}-\frac{l^{2}}{r^{3}}+\mathcal{O}\left(\frac{1}{r^{4}}\right).$
(63)

Near the horizon, we denote the small positive parameter $\epsilon=r-r_{h}$.
Then expanding the vector near the horizon, we obtain

$\displaystyle\phi^{r}(\theta=\frac{\pi}{2},\epsilon\rightarrow
0)=\frac{\sqrt{(M^{2}-a^{2})(l^{2}+4M^{2})}\left(M-\sqrt{M^{2}-a^{2}}\right)}{4a^{2}M^{2}}+\mathcal{O}(\epsilon^{\frac{1}{2}}).$
(64)

Note that the zero order is positive. Since
$\phi^{\theta}(\theta=\frac{\pi}{2})=0$, the direction of the vector is right
both at the horizon and large $r$. In particular, taking $a=0$, we obtain

$\displaystyle\phi^{r}(\theta=\frac{\pi}{2},r\rightarrow
r_{+})=\frac{1}{4M}\sqrt{\frac{l^{2}}{4M^{2}}+1}.$ (65)

For the Schwarzschild black hole, one has $\kappa=1/4M$ and
$g_{\varphi\varphi}^{H}=4M^{2}$, so this exactly produces the result given in
(32).

In summary, we conclude that the asymptotic behaviors of the vector for the
Kerr black hole are the same as that of the generic discussion given above.
And based on this result, we must have the fact that total topological number
of the TCOs in the Kerr black hole background is zero.

###  IV.3 Topological charge and stability

Here we would like to calculate the topological charge for the TCOs in the
Kerr black hole background. In order to show the characteristic behavior of
the unit vector $n$, we take the angular momentum $l$=1.5, 1.6827, and 1.9,
respectively. For the negative value of $l$, the behaviors are quite similar.
For simplicity, we take $a=0.98$ and $M=1$.

For $l$=1.5, we exhibit the unit vector $n$ in Fig. 5, where the arrows denote
the direction of the vector. It is clear that for different $r$ and $\theta$,
the direction of the vector is different. In the equatorial plane, one can see
that the vector is towards to right, and near $\theta$=0 and $\pi/2$, the
vector is up and down as expected. More importantly, there is no zero point
for the vector. Thus, the topological number defining as the sum of the
winding number for all the zero points is zero, i.e., $W=0$.

Figure 5: The blue arrows represent the direction of the unit vector $n$ for
$M$=1, $a=0.98$, and $l=1.5$.

If we set $l=1.6827$, the MSCO or ISCO will be present. For this particular
case, we show the unit vector $n$ in Fig. 6. Its behavior is similar to that
of Fig. 5. However, one IP1 point can be discovered at $r=1.61$, which
actually relates to the MSCO of the black hole. To calculate the winding
number for the MSCO, we construct a closed loop parametrized by the following
parametrized form

$\displaystyle\left\\{\begin{aligned} r&=c_{1}\cos\psi+c_{0},\\\
\theta&=c_{2}\sin\psi+\frac{\pi}{2}.\end{aligned}\right.$ (66)

For the closed loop $C_{1}$, we have $(c_{0},c_{1},c_{2})$=(1.61, 0.2, 0.6).
Along the loop $C_{1}$, we calculate the deflection angle $\Omega(\psi)$
defined as Wei2019

$\Omega(\psi)=\int_{C_{1}}\epsilon_{ab}n^{a}\partial_{i}n^{b}dx^{i}.$ (67)

The result is listed in Fig. 6. With the increase of $\psi$ from 0 to $2\pi$,
$\Omega(\psi)$ first increases, then decreases, and finally increases.
However, it does not make one complete loop. So the winding number for the
MSCO is zero, and thus $W=0$.

Figure 6: (a) The blue arrows represent the direction of the unit vector $n$
for $M$=1, $a=0.98$, and $l=1.6827$. (b) The red arrows represent the
direction of unit vector $n$ for $M$=1, $a=0.98$, and $l=1.9$. (c) The
deflection angle $\Omega(\psi)$ along the closed loop $C_{1}$. (d) The
deflection angle $\Omega(\psi)$ along the closed loops $C_{2}$ (red solid
curve) and $C_{3}$ (blue dashed curve). “IP” and “ZP” denote the locations of
the ISCO and the zero points of the vector.

As a third characteristic example, we take the angular momentum $l=1.9$. The
unit vector $n$ is exhibited in Fig. 6. Quite different from previous cases,
there are two zero points ZP2 and ZP3 located at $r=$1.36 and 2.46. In order
to calculate their winding numbers, we construct two closed loops C2 and C3
with parametrized coefficients $(c_{0},c_{1},c_{2})$=(1.36, 0.12, 0.6) and
(2.46, 0.4, 0.2), respectively. Along these loops, we show the deflection
angle in Fig. 6. With the increase of $\psi$, $\Omega(\psi)$ decreases for C2
and increases for C3. When $\psi$ vary from 0 to $2\pi$ denoting a complete
loop, we find that the vector makes one loop clockwise along C2, and
counterclockwise along C3. Therefore, we obtain $w_{ZP_{2}}$=-1 and
$w_{ZP_{3}}$=1 for the zero points ZP2 and ZP3. Considering that the TCOs with
small and large radii, respectively, are unstable and stable, one can conclude
that the unstable TCOs have negative winding number, while the stable ones
have positive winding number. Since the topological number is the sum of the
winding numbers of all the zero points, we still have $W=-1+1=0$. This
indicates that the total topological number still holds unchanged for this
case even when two zero points are present.

###  IV.4 Evolution of control parameter

In the above, we have shown several characteristic behaviors of the vector for
different values of the angular momentum acting as a control parameter. Here
we would like to consider the TCOs and their topological number under the
evolution of the control parameter.

We plot the radii of the TCOs as a function of the angular momentum $l$ in
Fig. 7, which can vary from $-\infty$ to $+\infty$. There are four different
branches. The signs “$+$” and “$-$” denote that the branches have positive and
negative winding numbers, respectively. For sufficiently negative angular
momentum, there are two branches of TCOs. The upper branch has a positive
winding number, while the lower one has a negative winding number; their sum
vanishes. With the increase of the angular momentum, these two branches
approach each other, and coincide at BP1 at $l$=-4.22. In the region
$l\in(-4.22,1.68)$, no TCO branch can be observed, and thus the topological
number is zero. Further increasing the angular momentum, we observe that two
branches extend from the point BP2 at $l$=1.68. The upper and the lower ones
have positive and negative winding numbers, respectively. The sum of the
winding number is still zero.

Figure 7: (a) The evolution of the radius $r_{\text{t}}$ of the TCO as the
angular momentum $l$. BP1 and BP2 are two bifurcation points. The signs “+”
and “-” denote that the branches have positive and negative winding numbers,
respectively. (b) The total topological number $W$ (black solid line) relating
to the TCOs as a function of the angular momentum $l$. Dashed lines denote the
winding numbers of these TCO branches. These black dots are for the
bifurcation points (BP). The black hole parameters are set to $M$=1 and
$a=0.98$.

From this behavior of $r_{\text{t}}$ of the TCOs, we can find that the
topological number always stays at zero for different values of the angular
momentum $l$. We show the details in Fig. 7. The dashed lines denote the
winding number of different TCO branches, and the black solid line is the
total topological number of the TCOs. It implies that the total topological
number always vanishes.

Moreover, from Fig. 7, one can clearly see that the points BP1 and BP2 are two
bifurcation points. Actually, they are the MSCOs for the Kerr black hole with
$a$=0.98 and $M$=1. Near these two points, we have

$\frac{d^{2}l}{dr_{\text{t}}^{2}}=-0.05\quad\text{and}\quad 2.31,$ (68)

which indicates that BP1 and BP2 are the annihilated and generated points,
respectively. This is also consistent with the result shown in Fig. 7.

Before ending this section, we briefly summarize our results. For the Kerr
black hole, the topological number relating to the TCOs vanishes, indicating
if the TCOs exist, they always come in pairs for the fixed angular momentum.
Note that although we here take $a$=0.98 and $M$=1 for an example, this result
holds for other values of the black hole mass and spin.

##  V Control parameter and topology

As shown above, we take the angular momentum $l$ as the control parameter to
construct the topology. Since the energy $E$ is also an alternative parameter
for the massive test particle, it is interesting to consider this case.

Here we need to reformulate the effective potential $\mathcal{V}$ as

$\displaystyle\mathcal{V}=-\frac{g_{tt}}{B}(l-l_{1})(l-l_{2}),$ (69)

where $l_{1,2}$ are given by

$\displaystyle l_{1}$ $\displaystyle=$
$\displaystyle\frac{-Eg_{t\phi}+\sqrt{B}\sqrt{g_{tt}+E^{2}}}{g_{tt}},$ (70)
$\displaystyle l_{2}$ $\displaystyle=$
$\displaystyle-\frac{Eg_{t\phi}+\sqrt{B}\sqrt{g_{tt}+E^{2}}}{g_{tt}}.$ (71)

It is obvious that $l_{1}\leq 0$ corresponds to the retrograde orbits.
Otherwise, $l_{2}$ corresponds to prograde orbits. So the conditions of the
TCO reduces to

$\displaystyle l=l_{1},\quad\text{and}\quad\partial_{r}l_{1}=0,$ (72)

for retrograde cases, and

$\displaystyle l=l_{2},\quad\text{and}\quad\partial_{r}l_{2}=0,$ (73)

for prograde cases. As a result, we need to construct two different vectors to
describe the TCOs for each rotating sense.

Case I: the vector $\phi_{1}$ is constructed as

$\displaystyle\phi_{1}^{r}=\frac{\partial_{r}l_{1}}{\sqrt{g_{rr}}},\quad\phi_{1}^{\theta}=\frac{\partial_{\theta}l_{1}}{\sqrt{g_{\theta\theta}}}.$
(74)

Via solving the zero point of $\phi_{1}$, we obtain that the corresponding
energy is $E=E_{-}$ given in (45). Plugging $E_{-}$ into $l_{1}$, one can
obtain the angular momentum $l=l_{-}$ as expected.

Case II: the vector $\phi_{2}$ reads

$\displaystyle\phi_{2}^{r}=\frac{\partial_{r}l_{2}}{\sqrt{g_{rr}}},\quad\phi_{2}^{\theta}=\frac{\partial_{\theta}l_{2}}{\sqrt{g_{\theta\theta}}}.$
(75)

For such case, we can obtain the zero point of $\phi_{2}$ at $E=E_{+}$ and
$l=l_{+}$.

Figure 8: Topological behaviors for the retrograde TCOs. (a) The evolution of
the radius $r_{\text{t}}$ of the TCO as the control parameter $E$. BP1 is a
bifurcation point at $E=0.96$. The signs “+” and “-” denote that the branches
have positive and negative winding numbers, respectively. (b) The total
topological number $W$ (black solid line) relating to the TCOs as a function
of the energy $E$. Dashed lines denote the winding numbers of these TCO
branches. The first and second black dots are for the bifurcation point and
phase transition point. The black hole parameters are set to $M$=1 and
$a=0.98$.

Figure 9: Topological behaviors for the prograde TCOs. (a) The evolution of
the radius $r_{\text{t}}$ of the TCO as the control parameter $E$. BP2 is a
bifurcation point at $E=0.77$. The signs “+” and “-” denote that the branches
have positive and negative winding numbers, respectively. (b) The total
topological number $W$ (black solid line) relating to the TCOs as a function
of the energy $E$. Dashed lines denote the winding numbers of these TCO
branches. The first and second black dots are for the bifurcation point and
phase transition point. The black hole parameters are set to $M$=1 and
$a=0.98$.

For each rotating sense, we can follow the treatment shown in Sec. III to
obtain the topological number by analyzing the asymptotic behaviors of the
vectors $\phi_{1}$ and $\phi_{2}$. Owing to that the process is quite similar,
we just present our result here. For both the retrograde and prograde TCOs, we
observe that, the topological number $W$ vanishes for $E\in(0,1)$, while when
the energy is larger than one, the topological number is $W$=-1. This is
mainly because that the direction of vectors $\phi_{1,2}^{r}$ at large $r$
changes when the energy crosses $E=1$. So there is a topological phase
transition at $E=1$ for both the rotating senses. As we know, in an
asymptotically-flat spacetime, the motion of the massive particles with energy
$E<1$ is bounded and cannot reach the infinity. Instead, these particles with
$E>1$ still keep moving at the infinity, and bounded orbits cannot be
observed. So the critical case is that these particles with energy $E=1$ can
stay at rest in the infinity. As a result, the energy $E$=1 plays a
characteristic point to distinguish these particles which cannot reach or
freely move at infinity. It also measures whether the bounded orbits exist or
not. From our topological approach, we clearly see that this physical feature
can be well-described by the topological phase transition at $E=1$.

In order to clearly show them, we exhibit our results in Figs. 8 and 9 for the
Kerr black hole with $M$=1 and $a$=0.98 for these two rotating senses. In
Figs. 8 and 9, we plot the evolutions of the radii $r_{\text{t}}$ of the TCOs
as the control parameter $E$ for the retrograde and prograde cases,
respectively. For the retrograde case, when $E<0.96$, there is no TCO branch,
and thus the total topological number is zero. When $E=0.96$, two branches of
negative and positive winding numbers emerge from the point BP1, and for the
prograde case, we observe a similar behavior near BP2 at $E=0.77$.

Near these two points BP1 and BP2, we have

$\frac{d^{2}E}{dr_{\text{t}}^{2}}=0.0018\quad\text{and}\quad 0.7631,$ (76)

which indicates that all these bifurcation points are generated points. Note
that the bifurcation points BP1 and BP2 are exactly the ones shown in Fig. 7.

As the energy approaches one, these two branches possessing positive winding
number tend to infinity $r_{\text{t}}\rightarrow\infty$. Nevertheless, the
total topological number always stays at zero as expected. When the energy is
larger than one, only one branch is left for each rotating sense, but both
them have negative winding number $w=-1$. Therefore, the topological number is
$W=-1$ for each case. These behaviors suggest that there is a topological
phase transition at $E=1$, where the topological number $W$ turns from 0 to
-1. Note that the values of the winding number correspond the stability of the
TCOs. As we have shown, the TCO with a negative winding number are unstable.
Any tiny perturbation will make the particle deviate from the TCO, so for the
particles with $E>1$, only unstable TCOs exist. Although these TCOs are very
limited for the study of astronomical observable effects, it can help us
understand the complete structure of TCOs in the black hole backgrounds.

We show the topological number $W$ in Figs. 8 and 9 for the retrograde and
prograde cases, respectively. It is clear that the topological number stays at
zero when $E<1$ for both cases, regardless of the bifurcation point. When
$E>1$, the topological number turns to -1 as expected. Thus, there is a
topological phase transition at $E$=1.

Noting that in our discussion above, we only take $a=0.98$. It is worth
examining whether topological pattern and phase transition depend on the mass
and spin of the black hole. As shown in Sec. IV.2, even though the black hole
mass and spin shall modify the coefficients of the expansion of the vector,
they do not affect the asymptotic behaviors. Thus, the total topological
number is independent of the black hole mass and spin. Nevertheless, the
locations of the TCOs and the MSCOs will be shifted. On the other hand, the
topological phase transition at $E$=1, is only dependent of the behavior of
$\phi^{r}$ at the infinity, so such a phase transition holds for arbitrary
values of the black hole mass and spin.

##  VI Conclusions and discussions

In this paper, we constructed a well-behaved topology for the TCOs in a
generic black hole spacetime, even though the TCOs are closely dependent of
the energy and angular momentum of the massive particle.

At first, we treated the angular momentum as the control parameter. Similar to
the treatment of the photon Cunhaa ; Cunhab , we started with the timelike
geodesics and constructed a vector $\phi$ by using the effective potential.
The TCO exactly locates at the zero point of $\phi$. Thus, one can endow each
of them with a local topological charge-the winding number, for the TCO.
Moreover, we studied the asymptotic behaviors of $\phi$ near the boundaries in
the $r$-$\theta$ plane. Via examining the axis limit, horizon limit, and
asymptotic limit, we found the value of the angular momentum does not change
the asymptotic behaviors of the vector. This strongly implies that the
topological approach can be applied to the TCOs, just like that for the LRs.

Then we considered the global and local topologies relating with the TCOs.
From the local viewpoint, we found that the positive or negative winding
number, respectively, corresponds to stable or unstable TCO. In particular,
the MSCO is just the bifurcation point that the stable and unstable TCOs can
generate from or annihilate into. Nevertheless, the sum of the winding number
stays the same before or after the bifurcation point. From the global
viewpoint, we observed that the total topological number, defined as the sum
of the winding numbers of all the zero points always vanishes for a generic
black hole spacetime. It suggests that the TCOs with the fixed angular
momentum come in pairs with one being stable and another being unstable. It is
worth pointing out that this result is universal and is independent of the
angular momentum of the particle and the black hole parameters. However, the
total topological number for the LR is -1, indicating there exists a standard
LR at least for a generic black hole.

As a specific example, we studied the Kerr black hole. As expected, the vector
is outwards at $\theta$=0 and $\pi$. Also the vector is to the right both at
$r=r_{\text{h}}$ and $\infty$. This pattern indicates that the total
topological number of the TCO is $W$=0 for the Kerr black hole. So, if one
treats the angular momentum as the control parameter, the stable and unstable
TCOs always come in pairs. Of particular interest is that the two MSCOs for
both the rotating senses act as the bifurcation points. The one with negative
angular momentum is an annihilated point and the other one with positive
angular momentum is a generated point. Although for different values of the
angular momentum, there may two or zero TCO branches, the total topological
number is always zero as expected.

Besides the angular momentum, the energy of the test particle is also an
alternative control parameter. Hence, we examine the corresponding topology
for both the rotating senses, respectively. The results show that the total
topological number is $W$=0 for $0<E<1$, while $W$=-1 for $E>1$. This suggests
that there exists a topological phase transition at $E$=1 both for the
retrograde and prograde cases, which actually measures different motions of
the massive particles at infinity. Moreover, the MSCOs in this case exactly
correspond to two generated points.

Before ending this paper, we would like to address three points for the
topology of the TCO that are different from the LR: i) For the LR, the energy
and angular momentum are linear, while for the TCO, they have a nonlinear
relation; ii) When the black hole background is fixed, the radius of the LR is
independent of the parameters of the photon, while for the TCO, its radius
closely depends on the test particle. It appears that the degeneracy of the
circular orbit is broken when the test particle gets mass, see Fig. 10 for an
example; iii) The total topological number is $W$=-1 for LR, while 0 for TCO.
This implies that for a generic black hole, there exists at least one standard
LR, while the TCOs always come in pairs for the fixed angular momentum. Here
we developed a general topological treatment for the TCOs, and applied it to
the Kerr black hole. We believe that more information about the TCOs could be
disclosed for other rotating and nonrotating black holes in GR or modified
theories of gravity by employing this topological approach, and more general
black hole solutions without $\mathcal{Z}_{2}$ will be further investigated.

## APPENDIX: REMOVAL OF DEGENERACY OF THE CIRCULAR ORBIT

Figure 10: Schematic illustration of the removal of degeneracy of the circular
orbit with $M$=1 and $a$=0.98. The PS of the Schwarzschild black hole is
broken by the black hole spin $a$ and particle mass $\mu^{2}$. Light red and
blue color regions are for the stable and unstable TCOs.
($r_{ISCO}^{a^{\pm}}$, $l_{ISCO}^{a^{\pm}}$)=(1.61, 1.68) and (8.94, -4.22).

Taking $M$=1 and $a=0.98$ as an example, we shown the removal of degeneracy of
the circular orbit due to the black hole spin and mass of the test particle in
Fig. 10. For the Schwarzschild black hole, the radius of the photon sphere
$r_{PS}$=3. When the black hole spin is present it is split into two; the
prograde and retrograde LRs located at $r_{LR}$=1.24 and 3.98, respectively.
Here $a^{\pm}$ denotes that the black hole spin $a$ and angular momentum $l$
have the same or opposite direction. When the test particle gains mass, its
motion will be described by the timelike geodesics rather than the null
geodesics. Each LR will turn into a continuous set of TCOs with the radius
extending from the value of LR to infinity. The unstable and stable TCOs are
separated by the MSCOs or ISCOs. Note that the radii of the unstable TCOs are
bounded by the ISCOs and LRs. For the $a^{+}$ branch, the particle angular
momentum starts at the MSCO with $l$=1.68 to $r_{TCO}^{a^{\pm}}$=1.24 and
$\infty$ with $l=\infty$, and for the $a^{-}$ branch, the particle angular
momentum starts at the MSCO with $l$=-4.22 to $r_{TCO}^{a^{\pm}}$=3.98 and
$\infty$ with $l=-\infty$.

## Acknowledgements

The authors would like to thank Dr. Jian-Zhi Yang for the useful discussions.
This work was supported by the National Natural Science Foundation of China
(Grants No. 12075103, No. 11875151 and No. 12247101).

## References

  * (1) P. V. P. Cunha, E. Berti, and C. A. R. Herdeiro, Light ring stability in ultra-compact objects, Phys. Rev. Lett. 119, 251102 (2017), [arXiv:1708.04211 [gr-qc]].
  * (2) P. V. P. Cunha and C. A. R. Herdeiro, Stationary Black Holes and Light Rings, Phys. Rev. Lett. 124, 181101 (2020), [arXiv:2003.06445 [gr-qc]].
  * (3) B. P. Abbott _et al_. (Virgo, LIGO Scientific), Observation of Gravitational Waves from a Binary Black Hole Merger, Phys. Rev. Lett. 116, 061102 (2016), [arXiv:1602.03837 [gr-qc]].
  * (4) K. Akiyama _et al_. (Event Horizon Telescope), First M87 Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole, Astrophys. J. 875, L1 (2019), [arXiv:1906.11238 [astro-ph.GA]].
  * (5) K. Akiyama _et al_. (Event Horizon Telescope), First Sagittarius A* Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole in the Center of the Milky Way, Astrophys. J. Lett. 930, L12 (2022).
  * (6) V. Cardoso, L. C. B. Crispino, C. F. B. Macedo, H. Okawa, and P. Pani, Light rings as observational evidence for event horizons: long-lived modes, ergoregions and nonlinear instabilities of ultracompact objects, Phys. Rev. D 90, 044069 (2014), [arXiv:1406.5510 [gr-qc]].
  * (7) S. Hod, On the number of light rings in curved spacetimes of ultra-compact objects, Phys. Lett. B 776, 1 (2018), [arXiv:1710.00836 [gr-qc]].
  * (8) M. Guo and S. Gao, Universal Properties of Light Rings for Stationary Axisymmetric Spacetimes, Phys. Rev. D 103, 104031 (2021), [arXiv:2011.02211 [gr-qc]].
  * (9) H. C. D. L. Junior, P. V. P. Cunha, C. A. R. Herdeiro, and L. C. B. Crispino, Shadows and lensing of black holes immersed in strong magnetic fields, Phys. Rev. D 104, 044018 (2021), [arXiv:2104.09577 [gr-qc]].
  * (10) R. Ghosh and S. Sarkar, Light rings of stationary spacetimes, Phys. Rev. D 104, 044019 (2021), [arXiv:2107.07370 [gr-qc]].
  * (11) H. C. D. L. Junior, J.-Z. Yang, L. C. B. Crispino, P. V. P. Cunha, and C. A. R. Herdeiro, Einstein-Maxwell-dilaton neutral black holes in strong magnetic fields: topological charge, shadows and lensing, Phys. Rev. D 105, 064070 (2022), [arXiv:2112.10802 [gr-qc]].
  * (12) Y. S. Duan, The structure of the topological current, Stanford Linear Accelerator Center (SLAC), Report No. SLACPUB- 3301, (1984).
  * (13) S.-W. Wei, Topological charge and black hole photon spheres, Phys. Rev. D 102, 064039 (2020), [arXiv:2006.02112 [gr-qc]].
  * (14) H.-S. Liu, Z.-F. Mai, Y.-Z. Li, and H. Lu, Quasi-topological Electromagnetism: Dark Energy, Dyonic Black Holes, Stable Photon Spheres and Hidden Electromagnetic Duality, Sci. China Phys. Mech. Astron. 63, 240411 (2020), [arXiv:1907.10876 [hep-th]].
  * (15) S.-W. Wei and Y.-X. Liu, Topology of black hole thermodynamics, Phys. Rev. D 105, 104003 (2022), [arXiv:2112.01706 [gr-qc]].
  * (16) P. K. Yerra and C. Bhamidipati, Topology of black hole thermodynamics in Gauss-Bonnet gravity, Phys. Rev. D 105, 104053 (2022), [arXiv:2202.10288 [gr-qc]].
  * (17) S.-W. Wei, Y.-X. Liu, and R. B. Mann, Black hole solutions as topological thermodynamic defects, Phys. Rev. Lett. 129, 191101 (2022), [arXiv:2208.01932 [gr-qc]].
  * (18) P. C. Peters, Gravitational radiation and the motion of two point masses, Phys. Rev. 136, B1224 (1964).
  * (19) V. Cardoso, C. F. B. Macedo, and R. Vicente, Eccentricity evolution of compact binaries and applications to gravitational-wave physics, Phys. Rev. D, 103, 023015 (2021), [arXiv:2010.15151 [gr-qc]].
  * (20) J. F. M. Delgado, C. A. R. Herdeiro, and E. Radu, Equatorial timelike circular orbits around generic ultracompact objects, Phys. Rev. D 105, 064026 (2022), [arXiv:2107.03404 [gr-qc]].
  * (21) A. Lehébel and V. Cardoso, The fate of observers in circular motion, Phys. Rev. D 105, 064014 (2022), [arXiv:2202.08850 [gr-qc]].
  * (22) 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), [arXiv:1609.01340 [gr-qc]].
  * (23) 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), [arXiv:1705.05461 [gr-qc]].
  * (24) J. Carot, Some developments on axial symmetry, Class. Quant. Grav 17, 2675 (2000).
  * (25) A. J. M. Medved, D. Martin, and M. Visser, Dirty black holes: Symmetries at stationary nonstatic horizons, Phys. Rev. D 70, 024009 (2004). [arXiv:gr-qc/0403026].
  * (26) C. A. R. Herdeiro and E. Radu, Kerr black holes with scalar hair, Phys. Rev. Lett. 112, 221101 (2014), [arXiv:1403.2757 [gr-qc]].
  * (27) J. F. M. Delgado, C. A. R. Herdeiro, and E. Radu, Rotating Axion Boson Stars, JCAP 06, 037 (2020), [arXiv:2005.05982 [gr-qc]].
  * (28) L.-B. Fu, Y.-S. Duan, and H. Zhang, Evolution of the Chern-Simons vortices, Phys. Rev. D 61, 045004 (2000), [arXiv:hep-th/0112033].


070004 2015 L. A. Pugnaloni G. Lumay, University of Liège, Belgium. 070004

This paper describes methods to measure flow rate and jamming onset of
granules discharged through a flat-bottom cylindrical hopper with a circular
orifice. The intrinsic jamming onset for ideal particles (spherical,
monodisperse, smooth) is experimentally measured by two independent methods,
with good agreement. For non-ideal particles, the normalized jamming onset
increases with elongated granule shape, broadened size distribution and
increased friction as measured by the drained angle of repose. An empirical
model of the jamming onset is introduced to quantify these effects over the
range of materials investigated. The jamming onset can be used as a measure of
differentiation between relatively free-flowing granules.

# Characterizing flowability of granular materials by onset of jamming in
orifice flows

Paul Mort[inst1] Email<EMAIL_ADDRESS>
(17 November 2014; 2 March 2015)

††volume: 7

99 inst1 Procter & Gamble Co., 5280 Vine Street, Cincinnati, OH 45217, USA.

##  1 Background

The motivation for this work is to find a means of characterizing and
differentiating flow quality of relatively free-flowing granular materials. On
the one hand, industry requires methods that can be performed relatively
simply and reproducibly. On the other hand, both industry and academia seek
better fundamental understanding of granular rheology via physical mechanisms
governing flows; in the case of industry, fundamental understanding should
extend to commercially-relevant granular materials. This paper correlates the
onset of jamming with granular characteristics, including particle size
distribution and particle shape. Details are provided on size and shape
characterization of commercially-relevant materials.

Characterization of cohesive powder flow has been relatively well established
using shear cells to quantify incipient bulk flow (i.e., yield loci) of
powders and granules. Flow functions calculated using yield loci can
differentiate between strongly-cohesive, mildly-cohesive and free-flowing
materials [1]; this is highly relevant to bulk storage and handling of
powders. However, shear cell measurements are relatively insensitive in
regards to differences among freely-flowing granules. For applications where
the quality of granular dynamic flow is of interest, we need a more sensitive
methodology, hence the motivation of the current work. Various approaches are
discussed in the literature, including using rotating drums [2, 3] and
impeller-driven flows [4, 5] to measure granular rheology. The current work
focuses on the onset of jamming in orifice flows to differentiate the quality
of otherwise freely-flowing granules.

Let us consider jamming in the context of a shear cell analysis. According to
continuum modeling, a “cohesionless” material measured using a shear cell
(i.e., a yield locus passing through the origin) implies flow through a hopper
cone opening of infinitesimal diameter. At this limit, the continuum model
fails because the orifice must be at least as large as the particle diameter.
Further, local packing and jamming effects require the opening to be
significantly larger than a single particle [6]. A theoretical analysis of
intrinsic cohesion suggests a boundary layer effect of 3-5 particle diameters
[7]; as a starting point, this can be taken as a theoretical minimum
dimensionless orifice size for an otherwise free-flowing granular material.
Experimental studies on the physics of jamming in granular flows report
jamming transitions on the basis of a dimensionless hopper opening at the
onset of the jamming transition [8, 9], where the dimensionless hopper opening
($D_{o}/d$) is defined by the ratio of the opening size ($D_{o}$) relative to
a monodisperse particle size ($d$). The current work is an empirical
investigation, seeking to elucidate the effects of both size distribution and
shape factors on the jamming probability of free-flowing granules.

##  2 Experimental

Figure 1: Schematic of a flat-bottom cylindrical hopper device (Flodex™).
Figure 2: Interpolation of Jamming Onset ($J$) at 25% mass discharged.

The experimental device used in this work was a flat-bottom cylindrical hopper
with a set of interchangeable disks having a range of discrete orifice sizes;
it is available as a commercial flow testing system, Flodex™(Hanson Research,
Chatsworth, CA, USA), with additional disks machined as required to extend the
range of orifice diameters. The range of orifice sizes used with the Flodex
included: 2.0, 2.5, 3.0, 3.5, 4, 5, 6, 7, 8, 9, 10, 12, 14, 16, 18 and 20 mm.
Experiments were done with common materials including narrowly-classified
glass beads, Ottawa sand, and a variety of granular detergent samples. In all
cases, these samples were relatively free-flowing according to shear cell
measurements.

A schematic of the orifice flow instrument is shown in Fig. 1. Approximately
100 ml of material is used in the test, filling the 5.7 cm diameter stainless-
steel hopper to a height of about 4 cm by pouring the granular sample through
the loading funnel. After the sample settles, the spring-loaded discharge gate
is opened and the sample is allowed to drain through the orifice into a
receiving cup below. Once the flow stops and remains stopped for 30 seconds,
the mass of discharged material is measured. Clogging is defined as a
persistent jammed state where the orifice remains obscured by the granules at
the point of flow stoppage. For each measurement, the mass% discharged is
calculated according to the formula: (mass% discharged) = 100 $\times$ (mass
discharged) / (sample mass). The average of the three mass% discharge
calculations is plotted as a function of the dimensionless orifice size
($D_{o}/d$), with the mass% discharged on the ordinate and the dimensionless
orifice size on the abscissa. This procedure is repeated using incrementally
larger orifice sizes until the hopper discharges without clogging for three
consecutive trials. The averaged data are linearly interpolated to find the
Jamming Onset ($J$), defined here as the value of the dimensionless orifice
size at the point of 25 mass% average discharge (Fig. 2).

For trials that drain without clogging, the drained angle of repose,
$\phi_{d}$, is calculated according to Eq. (1), derived assuming the granular
material discharges in the form of an inverted cone with a base diameter equal
to the hopper diameter, $D_{h}$, and truncated where the cone’s apex protrudes
through the orifice with diameter, $D_{o}$. Assuming the remaining volume of
retained material has cylindrical symmetry (Fig. 3), simple solid geometry and
the material’s repour (loose) bulk density, $\rho_{bulk}$, are sufficient to
relate the measured mass of retained material, $M_{ret}$, to the drained angle
of repose.

$\phi_{d}=\arctan\left[\frac{24M_{ret}}{\pi\rho_{bulk}(2D_{h}^{3}-3D_{o}D_{h}^{2}+D_{o}^{3})}\right]$
(1)

Orifice size is made dimensionless by scaling the orifice size ($D_{o}$) to a
characteristic particle size ($d$). The Sauter mean ($d32$) is used as the
characteristic size for this study. Note the Sauter mean is weighted toward
the surface area of the granules, the surface area being relevant to inter-
particle contact and frictional interactions of the flow. Studies of cohesive
dry powders show correlation between the Sauter mean size and flow properties
[10]. Size distributions were measured by sieving, fitting the size
distribution using a log-normal distribution model to obtain a mass-based
geometric mean ($d43$) and geometric standard deviation ($\sigma_{g}$) for
each sample, then converting to the Sauter mean using the Hatch-Choate
relation [11], Eq. (2). Examples of size distribution characteristics are
shown in Fig. 4, for select samples A-D discussed in detail later in the
results section.

Figure 3: Vertical cross section of residual mass in a cylindrical hopper
after draining successfully without a clog.  Figure 4: Examples of particle
size distribution analysis for four samples, A-D. The data points are from
sieve analysis; the fitted curves are based on a log-normal distribution,
weighted to the mass on each sieve.  Figure 5: Image analysis freeze frame and
cross-sectional aspect ratio (AR) analysis: a) grayscale image; b) B&W image
resulting from threshold analysis with particles outlined in red, $AR=d1/d2$.
Figure 6: Discharge rate data for Sample D, measured over a series of orifice
diameters (4 mm to 9 mm). For each condition, the flow rate is the time
derivative of mass discharged. The flow rate data are plotted according to the
Beverloo relation in Fig. 9.

$\ln{d32}=\ln{d43}-\ln^{2}{\sigma_{g}}$ (2)

Shape characteristics were measured using an automated image analysis system,
Solids Sizer (JM Canty, Buffalo, NY), using statistical averaging of $10^{4}$
cross-sectional particle images per sample. An example of particles counted in
a freeze-frame captured by the image analysis system is given in Fig. 5. Note
that the particles are in free-fall at the point of image capture, so the
cross sectional image is randomly distributed over the possible states of
rotation. The median aspect ratio ($AR50$) was used as the characteristic
shape factor in this study.

The granular discharge rate was measured by placing the receiver cup on an
electronic balance with data acquisition to record the discharge mass with
time. This experiment was done using orifice sizes above the jamming onset.
Discharge data are provided in Fig. 6, using Sample D as an example. Sample D
is discussed in more detail in the results section. The flow rate was
determined by taking the time derivative (slope) of the cumulative discharge.
The mass flow rate is converted to a volumetric rate using the repour bulk
density.

##  3 Results and Discussion

###  3.1 Analysis of jamming onset using a multivariate model

Results for nearly “ideal” materials including tightly classified glass beads
and washed Ottawa sand are shown on Fig. 7. These samples are nearly mono-
disperse. The glass beads are nearly spherical, while the sand is slightly
irregular in shape. Both have relatively smooth-surfaces and a low drained
friction angle. The jamming onset for the glass beads is in the range of 3 to
5 particle diameters. The washed sand has a slightly higher jamming onset of
about 7 diameters. Other samples characterized in this work included
commercial granules having similar compositions, but with a range of particle
size and shape characteristics, the effect of which resulted in relative
jamming onsets up to 30 and drained friction angles in excess of 60 Degrees.
Because the commercial samples have several aspects of variability, a
multivariate model was formulated to assess the relative effects of friction,
size distribution and shape characteristics.

Figure 7: Cylindrical hopper discharge data, jamming onset calculations and
drained angle of repose for glass bead and washed sand samples. Figure 8:
Multivariate regression of jamming onset model [Eq. (3)]; deviation from the
diagonal represents uncertainty in the empirical model. Solid data points
labeled A-D indicate samples used in Flow Rate analysis (Fig. 9); their size
distributions are shown in Fig. 4.  Figure 9: Regression [Eq. (5)] of flow
rate data collected for selected samples (A-D) having a range of jamming
onsets.

Data from the broader set of 23 samples (Table 1), including commercial
granular materials and manipulations thereof, were used to generate a
multivariate power-law model for the jamming onset ($J$) as a function of
particle parameters [Eq. (3)]. The model parameters include the shape factor
(median aspect ratio, $AR50$), size distribution breadth (geometric standard
deviation, $\sigma_{g}$) and an excess friction factor defined as
$XF=\max[1,\tan(\phi_{d})]$. The range in breadth of size distribution was
from about 1.0 to 2.2; the range of median aspect ratio was from about 1.0 to
1.4. The fitted coefficients ($a$, $b$, $c$) represent exponents in the power-
law form of the model. While this equation is purely empirical, its form is
logical in the sense that it reduces to an intrinsic jamming intercept
($k_{j}$) when the other terms are minimized. In other words, extrapolation to
$k_{j}$ represents the jamming onset of an idealized sample (mono-disperse,
spherical, smooth, low friction).

$\displaystyle\ln(J)=$ $\displaystyle\ln(k_{j})+a\ln(AR50)$
$\displaystyle+b\ln(\sigma_{g})+c\ln(XF)$ (3)

The results of the multivariate model are shown graphically in Fig. 8. While
the correlation coefficient ($\sim$0.86) indicates that the model is
reasonably predictive, the scatter suggests uncertainty, perhaps in the choice
of model parameters, their interaction and/or uncertainty in measurements. For
example, the model may be over-simplifying the parameter space, ignoring
potentially important factors such as more complex shape factors, the breadth
of the shape distribution, and interactions between size and shape
distributions.

The statistical analysis is shown in Table 2. The most statistically
significant factors are the intrinsic jamming intercept, $k_{j}$, excess
friction factor, $XF$, and the geometric size distribution, $\sigma_{g}$. The
value of $k_{j}$ is $\sim$2.9, at the lower end of the range predicted by Wier
[7]. While the particle shape factor, $AR50$, is somewhat less significant
statistically, its high coefficient (2.38) indicates that particle shape can
have a strong impact on jamming, even over the relatively narrow range tested.

Table 1: Experimental data for relative jamming onset ($J$) as a function of excess friction factor ($XF$), geometric size distribution ($\sigma_{g}$), and median aspect ratio ($AR50$). Sample | $J$ | $XF$ | $\sigma_{g}$ | $AR50$
---|---|---|---|---
A | 18.6 | 1.49 | 1.84 | 1.36
B | 16.2 | 1.01 | 1.76 | 1.35
C | 11.8 | 1.05 | 1.63 | 1.35
D | 7.9 | 1.00 | 1.37 | 1.26
E | 6.8 | 1.00 | 1.33 | 1.30
F | 7.1 | 1.00 | 1.34 | 1.24
G | 8.9 | 1.10 | 1.64 | 1.25
H | 9.6 | 1.08 | 1.48 | 1.32
I | 28.9 | 1.90 | 1.79 | 1.41
J | 15.1 | 1.00 | 2.14 | 1.26
K | 18.3 | 1.53 | 1.50 | 1.32
L | 28.2 | 1.24 | 1.79 | 1.34
M | 15.6 | 1.00 | 1.56 | 1.41
N | 19.1 | 1.83 | 1.56 | 1.41
O | 6.5 | 1.08 | 1.47 | 1.28
P | 7.3 | 1.00 | 1.39 | 1.27
Q | 9.1 | 1.05 | 1.41 | 1.26
R | 7.6 | 1.09 | 1.32 | 1.27
S | 17.6 | 1.73 | 1.31 | 1.26
T | 6.8 | 1.00 | 1.23 | 1.21
U | 7.1 | 1.00 | 1.10 | 1.20
V | 4.3 | 1.00 | 1.05 | 1.10
W | 3.5 | 1.00 | 1.05 | 1.15
Table 2: Parameter estimates and standard errors obtained from multi-variate regression of Eq. (3) using the data of Table 1. “Prob$>|t|$” is the probability that the true parameter value is zero, against the two-sided alternative that it is not; values less than about $0.05$ are typically regarded as highly significant. Term | Estimate | Std Error | P$>|t|$
---|---|---|---
$\ln(k_{j})$ | 1.072 | 0.225 | 0.0001
$a$ | 2.380 | 1.267 | 0.0758
$b$ | 1.404 | 0.398 | 0.0023
$c$ | 1.095 | 0.265 | 0.0006

###  3.2 Flow rate measurements and Beverloo analysis

Flow rate data were analyzed using the Beverloo equation [12]. It is shown
here in volumetric form [Eq. (4)], where $V^{\prime}$ is the volumetric feed
rate, $D_{o}$ is the orifice diameter, $d$ is the characteristic particle size
(Sauter Mean), $g$ is the gravitational constant and $C$ and $k_{f}$ are fit
parameters; $D_{o}-k_{f}d$ represents a reduced orifice size for active flow
caused by a boundary layer that scales with doparticle size. For the purpose
of regression analysis as a function of relative orifice size ($D_{o}/d$), the
Beverloo equation can be written in linearized form, Eq. (5), where $C$ is
solved using the regression slope and $k_{f}$ using the intercept.

$V^{\prime}=C\sqrt{g}\left(D_{o}-k_{f}d\right)^{5/2}$ (4)

$\frac{V^{\prime 2/5}}{d}=C^{2/5}g^{1/5}\left(\frac{D_{o}}{d}-k_{f}\right)$
(5)

Several flow rate data sets are shown in Fig. 9. These data sets (Samples A-D)
represent similar materials, but with different shape and size distribution
characteristics. Sample A is a commercial granule; samples B and C are
manipulated by classification; and sample D is classified then further rounded
in a layering process. The jamming onsets ($J$) of these samples span a
significant range from about 8 to 20. However, the Beverloo regressions are
remarkably similar. All lie on a common slope ($C\sim 0.57$) and extrapolate
to a similar zero-flow intercept ($k_{f}\sim 2.64$). In other words, the
granular flow behavior is consistent for orifice sizes above the jamming
onset, but subtle differences in the granular characteristics and frictional
properties have a significant effect on the onset of jamming. This suggests
that the onset of jamming is a useful measure to differentiate the flow
behavior of relatively free-flowing granular materials.

Consider $k_{f}$ to represent an “intrinsic jamming onset” at which the flow
rate goes to zero. Note that the extrapolated value for intrinsic jamming
(i.e., no flow) is only slightly lower than the power-law intercept of the
multivariate model ($k_{j}\sim 2.9$). Recall, the measurement of the jamming
onset, $J$, is based on the interpolated orifice size for an average 25%
sample discharge from the test hopper. Thus, it is logical that $k_{j}$, which
is based on multivariate analysis of $J$, may slightly exceed $k_{f}$, which
implies zero flow. On the other hand, the slope of the jamming probability
function tends to be quite steep in the range from zero to 25% discharge. So,
for the purpose of this paper, we can merge the interpretation of $k_{j}$ and
$k_{f}$ —both represent an intrinsic jamming limit. Indeed, even though each
one was measured independently using different methods, their values are
nearly identical.

An interpretation of this result is that the “intrinsic jamming onset” is
characteristic of orifice flow for an ideal particle (monodisperse, spherical,
no excess friction), and that differences between the intrinsic onset and the
actual measured jamming onset are due, in large part, to irregularities in the
size distribution, shape and frictional roughness of the granules.

##  4 Conclusions

The statistical analysis of jamming is a promising approach for
characterization of relatively free-flowing granular materials. In this work,
we evaluated the effect of dimensionless orifice size on the onset of jamming
and flow rate. A multivariate analysis of the jamming onset suggests an
intrinsic jamming probability of about 2.9 particle diameters for an idealized
sample; irregularities such as size distribution, shape or excess friction
support jamming across an increased number of particle diameters. While all
three parameters are significant in the regression analysis, the result
suggests that shape (aspect ratio) is an especially important contributor to
jamming. Analysis of the flow rate data using the Beverloo equation
extrapolates to a similar value for the intrinsic jamming onset ($\sim$2.6
particle diameters), providing good experimental agreement between two
independent methods.

The experimental methods described herein are relatively simple. Continued
work covering a wider range of materials can help to build stronger empirical
models for the effects of shape, size distribution and other frictional
properties on jamming probability. In addition, first principle modeling of
dense flows may help to elucidate the underlying mechanisms of jamming and the
real effects of non-ideal particle characteristics, and perhaps help to build
a more theoretical model.

Lastly, the reader should temper the results presented with the caveat that
all of the samples used herein are relatively free-flowing granular materials.
As materials become more cohesive and/or compressible, it is not clear that
the flow rate and jamming onset experiments will continue to converge on a
common intrinsic jamming onset.

###### Acknowledgements.

The original data for this paper were collected and analyzed as part of a
student internship at Procter & Gamble by Derek Geiger, University of
Michigan, with the assistance of Mark Wandstrat at P&G, and presented at the
2007 Annual Meeting of the American Institute of Chemical Engineers [13]. This
work at P&G followed an earlier collaboration with Professor R. P. Behringer
[14]; the author acknowledges Professor Behringer’s guidance and insight on
using fundamental granular physics to address problems of industrial
relevance.

## References

  * [1] J A H de Jong, A C Hoffmann, H J Finkers, Properly determine powder flowability to maximize plant output, Chem. Eng. Prog. 95, 25 (1999).
  * [2] G Lumay, F Boschini, K Traina, S Bontempi, J-C Remy, R Cloots, N Vandewalle, Measuring the flowing properties of powders and grains, Powder Technol. 224, 19 (2012).
  * [3] A Alexander, B Chaudhuri, A Faqih, F Muzzio, C Davies, S Tomassone, Avalanching flow of cohesive powders, Powder Technol. 164, 13 (2006).
  * [4] M Kheiripour Langroudi, S Turek, A Ouazzi, G I Tardos, An investigation of frictional and collisional powder flows using a unified constitutive equation, Powder Technol. 197, 91 (2010).
  * [5] M Leturia, M Benali, S Lagarde, I Ronga, K Saleh, Characterization of flow properties of cohesive powders: A comparative study of traditional and new testing methods, Powder Technol. 253, 406 (2014).
  * [6] A H Birks, M S A Bradley, R Farnish, The conversion of the analytical simple shear model for the Jenike failure locus into principle stress space and implication of the model for hopper design, In: Handbook of Conveying and Handling of Particulate Solids, Eds. A Levy, H Kalman, Pag. 95, Elsevier (2001).
  * [7] G J Weir, The intrinsic cohesion of granular materials, Powder Technol. 104, 29 (1999).
  * [8] I Zuriguel, A Garcimartín, D Maza, L A Pugnaloni, J M Pastor, Jamming during the discharge of granular matter from a silo, Phys. Rev. E 71, 051303 (2005).
  * [9] K To, P-Y Lai, H K Pak, Jamming of granular flow in a two-dimensional hopper, Phys. Rev. Lett. 86, 71 (2001).
  * [10] D Geldart, E C Abdullah, A Verlinden, Characterisation of dry powders, Powder Technol., 224, 19 (2012).
  * [11] T Hatch, S P Choate, Statistical description of the size properties of non-uniform particulate substances, J. Franklin Inst., 207, 369 (1929).
  * [12] W A Beverloo, H A Leniger, J Van de Velde, The flow of granular solids through orifices, Chem Eng. Sci., 15, 260 (1961).
  * [13] D Geiger, M Wandstrat, P Mort, Granular flow through an orifice – effect of granule size and shape distributions, In: Proceedings of the AIChE Annual Meeting, Salt Lake City, UT, USA (2007).

  * [14] P Mort, K McKenzie, J Wambaugh, R Behringer, Granular flow through an orifice – effect of granule size and shape distributions, In: Proceedings of the World Congress of Particle Technology 5, Orlando, FL, USA (2006).