e1e-mail<EMAIL_ADDRESS>e2e-mail<EMAIL_ADDRESS>e3e-mail<EMAIL_ADDRESS>
11institutetext: Research Centre for Theoretical Physics and Astrophysics,
Institute of Physics, Silesian University in Opava, Bezručovo nám.13, CZ-74601
Opava, Czech Republic. 22institutetext: Department of Mathematics, COMSATS
University Islamabad, Lahore Campus, 54000 Lahore, Pakistan.

#  Charged particle dynamics in parabolic magnetosphere around Schwarzschild
black hole

Martin Kološe1,addr1  Misbah Shahzadie2,addr2  Arman Tursunove3,addr1

(Received: date / Accepted: date)

###### Abstract

The study of charged particle dynamics in the combined gravitational and
magnetic field can provide important theoretical insight into astrophysical
processes around black holes. In this paper, we explore the charged particle
dynamics in parabolic magnetic field configuration around Schwarzschild black
hole, since the paraboloidal shapes of magnetic field lines around black holes
are well motivated by the numerical simulations and supported by observations
of relativistic jets. Analysing the stability of bounded orbits and using the
effective potential approach, we show the possibility of existence of stable
circular off-equatorial orbits around the symmetry axis. We also show the
influence of radiation reaction force on the dynamics of charged particles, in
particular on the chaoticity of the motion and Poincaré sections, oscillatory
frequencies, and emitted electromagnetic spectrum. Applied to Keplerian
accretion disks, we show that in parabolic magnetic field configuration, the
thin accretion configurations can be either destroyed or transformed into a
thick toroidal structure given the radiation reaction and electromagnetic-disk
interactions included. Calculating the Fourier spectra for radiating charged
particle trajectories, we find that the radiation reaction force does not
affect the main frequency peaks, however, it lowers the higher harmonics
making the spectrum more flat and diluted in high frequency range.

###### Keywords:

black hole, magnetic field, charged particle, chaotic dynamics

††journal: Eur. Phys. J. C

August 28, 2024

##  1 Introduction

Long-range forces provided by gravitational and electromagnetic (EM)
interactions are of crucial importance for the proper understanding of high-
energy processes around black holes (BHs). There are convincing observational
evidences that magnetic fields (MFs) must be present in the vicinity of BHs
Eatough-etal:2013:Natur: ; Daly:2019:APJ: . Orders of magnitude of MFs around
BHs may vary from a few Gs up to $10^{8}$ Gs and more, depending on the source
generating the field. For stellar-mass BHs observed in X-ray binaries, the
characteristic strength of MFs is of the order of $10^{8}$ Gs, while for
supermassive BHs (SMBHs), the characteristic strength is of the order of
$10^{1}-10^{4}$ Gs Eatough-etal:2013:Natur: ; Gold-etal:2017:ApJ: ;
EHT:2021:ApJ: . Since the energy densities of MFs of such orders are not
enough to make a sufficient contribution to the geometry of the background
spacetime, in realistic astrophysical situations the spacetime metric around a
BH can be fully described by the Kerr or Schwarzschild solution of the
Einstein field equations.

The weakness of the MF around BH, in the sense that the MF does not affect the
spacetime geometry, is compensated by the large value of the charge-to-mass
ratio for elementary particles, whose motion will be essentially affected by
MFs already of the order of few Gs. For a charged test particle with charge
$q$ and mass $m$ moving in the vicinity of a BH with mass $M$ immersed in MF
of the strength $B$, one can introduce a dimensionless “magnetic parameter”
$\cal{B}$ reflecting the ratio between the Lorentz force (LF) and
gravitational force. For a relativistic electron orbiting a BH at the distance
of the event horizon scale, one can estimate this ratio as follows

${\cal{B}}\sim
B\,\frac{q}{m}\,\frac{GM}{c^{4}}\approx\,10^{11}\,\,\frac{q}{e}\,\frac{m_{\rm
e}}{m}\,\frac{B}{10^{8}{\rm G}}\,\frac{M}{10M_{\odot}},$ (1)

that is very large due to large values of the specific charge $q/m$ for
elementary particles. Therefore, the influence of the EM force on the motion
of charged particles can be dominant in realistic situations.

Figure 1: Parabolic MF configuration around Schwarzschild BH. Quasi neutral
Keplerian accretion disk formed by particles on circular orbits close to the
equatorial plane is restricted in its inner edge by the equatorial innermost
stable circular orbit (inISCO). The dashed curve is the boundary between the
accretion disk corona and jet funnel region close to the $z$-axis.

Within the theory of general relativity, the gravitational influence of a
central compact object (BH) can be well described by the Kerr or Schwarzschild
spacetime Stu-etal:2020:Universe: , but unfortunately, we don’t have such an
elegant and well accepted solution for EM field influence – instead, there
exist various models for BH magnetosphere depending on the complexity of
included physical effects. As a zero approximation to the problem of BH
magnetosphere one can use the vacuum solutions of Maxwell equations in the
curved background as a starting point. In such approximation there is no
source $\rho=0$, or EM currents $J^{\alpha}=0$ included, and the solution can
be considered as an EM field of isolated BH in vacuum. However, electrovacuum
solutions have limited astrophysical relevance, as they do not include the
effects of material orbiting around BH in the form of a plasma. In the force-
free approximation plasma is taken into account, but MF energy density is
still dominating over matter energy density $B^{2}\gg\rho c^{2}$. In such
models, the matter is used to generate EM currents inside the BH
magnetosphere, but the matter is still locked in the MF. Plasma effect on the
force-free magnetosphere of BH has been studied in the well-known work of
Blandford & Znajek Bla-Zna:1977:MNRAS: , where also the EM mechanism of
rotational energy extraction from BH has been proposed. Several numerical
techniques have been also employed, but the exact shape and intensity of the
BH magnetosphere are not yet properly resolved, although a strong connection
to the accretion processes is evident Meier:TheEngineParadigm ; Galishnikova-
etal:2022:XXX: . One can also use the general relativistic
magnetohydrodynamics (GRMHD) or general-relativistic particle-in-cell (GRPIC)
numerical simulations to test EM field configurations and try to construct a
relevant BH magnetosphere model. A simple and still relevant BH magnetosphere
model inspired by GRMHD and GRPIC simulations is the parabolic MF
configuration Nak-etal:2018:APJ: ; Por-etal:2019:APJS: ; Kol-Jan:2020:RAG: ;
Cri-etal:2020:AA: .

MF effects on the astrophysical phenomena occurring in the vicinity of BHs
cannot be neglected and the aim of this article is to find close connections
between a simple realistic model for the BH magnetosphere with a radiative
processes in the combined gravitational and EM fields. Simple charged particle
model has limited validity in description of matter around BH – one should
assume low particle density, i.e., that the mean time of two particle
interaction is longer than orbital time around BH. One should also note that
the typical Larmor radius of a particle in BH magnetosphere is much smaller
than the Schwarzschild radius. Such low particle density assumption is correct
only for the jet and for disks of low-luminosity active galactic nuclei (AGN),
like Sgr A* or M87, while plasma inside thin disk around luminous AGN or XRBs
may be quite collisional due to thin accretion disk’s high density.

Our work can be seen as preliminary basic study of realistic BH magnetosphere,
which should have according to the recent observation of synchrotron radiation
and numerical modelling parabolic character Nak-etal:2018:APJ: . The article
is structured as follows. We explore charged particle motion in parabolic BH
magnetosphere model using effective potential approach in section 2, then we
examine the influence of radiation reaction (RR) damping force on particle
chaotic dynamics in section 3. We tested neutral Keplerian disk ionisation in
section 4, estimate observable frequencies in section 5, and provide some
astrophysical application in section 6.

Throughout the paper, we use the spacelike signature $(-,+,+,+)$. Greek
indices are taken to run from $0$ to $3$, while Latin indices run only over
the spatial coordinates ($r,\theta,\phi$). However, for expressions having
astrophysical relevance, we use the physical constants explicitly.

Figure 2: Effective potential $V_{\rm eff}$ for charged particle dynamics as
2D function of Cartesian $x,z$ coordinates. We can see the effective potential
minima at the off-equatorial planes for ${\cal{B}}<0$ (top left) and the in-
equatorial minima for ${\cal{B}}>0$ (top right). The lower row of figures
shows the sections of the effective potential at given energy and the
corresponding minima. Value of the parameter $w$ has been set to $w=3/4$, but
similar behavior is observed for different values of $w\in(0,1)$.

##  2 Charged particles in black hole magnetosphere

The line element describing the spacetime of non-rotating Schwarzschild BH
with mass $M$ is given as

$\mathrm{d}s^{2}=-f(r)\mathrm{d}t^{2}+f^{-1}(r)\mathrm{d}r^{2}+r^{2}(\mathrm{d}\theta^{2}+\sin^{2}\theta\mathrm{d}\phi^{2}),$
(2)

where the function $f(r)$ takes the form

$f(r)=1-\frac{2M}{r}.$ (3)

Hereafter, we put $M=1$, i.e., we use dimensionless radial coordinate $r$ (and
time coordinate $t$).

In our magnetosphere model, we assume the axial symmetry and absence of
electric field, hence the only non-zero component of EM four-vector potential
$A_{\mu}$ will be $A_{\phi}$, so that $A_{\mu}=(0,0,0,A_{\phi})$. We will
consider heuristic BH magnetosphere model given by the parabolic MF solution
motivated by the GRMHD plasma simulations McK-Nar:2007A:MNRAS: ; McK-
Nar:2007B:MNRAS: ; Tch-Nar-McK:2010:APJ: ; Cri-etal:2020:AA:

$A_{\phi}(r,\theta)=\frac{B}{2}\,r^{w}\,(1-|\cos\theta|),$ (4)

where $B\in(-\infty,\infty)$ is MF intensity and $w\in[0,1.25]$ is declination
of MF lines. Such parabolic field configuration should represent magnetic
field inside BH jet region (close to $z$-axis) and more complex configuration
is expected in accretion disk region (close to equatorial plane). One should
note that (4) is a simplified version of the Blandford-Znajek (BZ)
paraboloidal model Bla-Zna:1977:MNRAS:  and it is not a solution of a vacuum
Maxwell equation in Schwarzschild spacetime. The absolute value in $A_{\phi}$
definition is making the solution to split in equatorial plane and preserving
MF divergence, so MF lines go into the BH below the equatorial plane and
outgoing above the equatorial plane, see Fig. 1. Such split magnetic
configuration would require existence of conductive plasma current sheet in
equatorial plane, thus a material orbiting around and falling into BH is
supporting this model. In the case with $w=1$ one gets a limit of BZ
paraboloidal model and for $w=0$ the BZ split monopole solution Bla-
Zna:1977:MNRAS: . The value, which is mostly used for BH magnetosphere in jet
funnel is $w=3/4$ Nak-etal:2018:APJ: .

The orthonormal components of the MF ${\bf
B}=(B^{\widehat{r}},B^{\widehat{\theta}},B^{\widehat{\phi}})$, measured by
zero-angular momentum observer (ZAMO) can be determined by Tur-et-
al:2021:PhRvD:

$B^{\widehat{i}}=\frac{1}{2}\eta_{ijk}\sqrt{g^{jj}g^{kk}}F_{jk},$ (5)

where $\eta_{ijk}$ is the Levi-Civita tensor in three dimensional space, and
$F_{ij}=A_{j,i}-A_{i,j}$ is the antisymmetric tensor of the EM field. The non-
vanishing orthonormal components of the parabolic MF take the form

$\displaystyle B^{\widehat{r}}$ $\displaystyle=$
$\displaystyle\frac{B}{2}\leavevmode\nobreak\
\frac{\cos\theta}{|\cos\theta|}\leavevmode\nobreak\ r^{w-2},$ (6)
$\displaystyle B^{\widehat{\theta}}$ $\displaystyle=$
$\displaystyle-\frac{B}{2}\left(\frac{1-|\cos\theta|}{\sin\theta}\right)\sqrt{f(r)}\leavevmode\nobreak\
r^{w-2}\leavevmode\nobreak\ w.$ (7)

The magnitude of MF B as a function of $r$ and $\theta$ takes the form

$\displaystyle|\textbf{B}|$ $\displaystyle=$
$\displaystyle\frac{1}{2}[B^{2}r^{2w-5}\,\\{r-\frac{1}{2}\left(4\leavevmode\nobreak\
|\cos\theta|-3-\cos 2\theta\right)\csc^{2}\theta$ (8) $\displaystyle\times$
$\displaystyle(r-2)\,w^{2}\\}]^{1/2}.$

The realistic accretion states from GRMHD simulations that get magnetically
saturated do form such a split monopole state that lives transiently until
accretion destroys it again (i.e. a thin current sheet at the base of a
paraboloidal jet that looks like a split monopole nearby the BH) Ripperda-
etal:2022:ApJ: . The realistic or simulated BH magnetospheres Tch-Nar-
McK:2010:APJ: ; Cri-etal:2020:AA:  have a more complicated structure. Despite
this, our simplified model (4) captures their most important features. One can
distinguish three different regions around an accreting BH (see Fig. 1): the
accretion disk/torus, where the matter density is high $\rho\sim\rho_{\rm
max}$ and matter dominates over the turbulent MF; corona region, where the
matter density is much lower $\rho\leq\rho_{\rm max}$, but the matter still
dominates over the still turbulent MF; jet funnel, where the matter component
is very weak $\rho\leq 10^{-6}\rho_{\rm max}$ and regular MF with parabolic
shape dominates the region. The MF lines arising from the BH event horizon in
equatorial plane are giving the boundary between the corona and jet funnel
region. Corona-jet boundary condition is independent of magnetic intensity
parameter $B$ and reads

$2^{w}-r^{w}\left(1-|\cos\theta|\right)=0.$ (9)

In our simplified Schwarzschild BH magnetosphere model (4), the MF lines are
not rotating and light cylinder will not be formed. Mostly due to the density
differences, the mean free path for charged particles differs dramatically in
all three (disk/corona/jet) regions. It is interesting to note that at least
for M87 and Sgr A* parameters, also the disk is collisionless due to the small
mass accretion rate and low plasma mass density. In this article, we will
assume the description of collisionless charged test particle dynamics can be
well applied at least in jet funnel and corona regions and hence we will focus
on these two regions leaving the description of accretion disk to a more
advanced GRMHD Tch-Nar-McK:2010:APJ: ; Nak-etal:2018:APJ:  or PIC Cri-
etal:2020:AA: ; Hir-etal:2021:APJ:  simulations.

Figure 3: For attracting LF (${\cal{B}}<0$), we plot the radial positions of
stable circular orbits for charged particles in the off-equatorial plane.
Different MF configurations $w=0,0.75,1$ have been plotted. The grey disk
shows the BH horizon, the numbers on each solid curve indicate the values of
the MF parameter ${\cal{B}}$, while the points of these curves show the
position of stable off-equatorial circular orbits. The dotted curve represents
the offISCO for charged particles in the off-equatorial plane, while the
dashed curve shows the corona/jet boundary.

The charged test particle motion is described by the covariant Lorentz
equation

$m\frac{Du^{\mu}}{\mathrm{d}\tau}=qF^{\mu}_{\nu}u^{\mu},$ (10)

where $u^{\mu}$ is the four-velocity of the particle with the mass $m$ and
charge $q$, normalized by the condition $u^{\mu}u_{\mu}=-1$, $\tau$ is the
proper time of the particle, and $F_{\mu\nu}=A_{\nu,\mu}-A_{\mu,\nu}$ is the
antisymmetric tensor of the EM field. Using Hamiltonian formalism for the
charged particle motion, we can write

$H=\frac{1}{2}g^{\alpha\beta}(\pi_{\alpha}-qA_{\alpha})(\pi_{\beta}-qA_{\beta})+\frac{1}{2}\,m^{2},$
(11)

where the kinematical four-momentum $p^{\mu}=mu^{\mu}$ is related to the
generalized (canonical) four-momentum $\pi^{\mu}$ by the relation

$\pi^{\mu}=p^{\mu}+qA^{\mu},$ (12)

that satisfies the Hamilton equations in the form

$\frac{\mathrm{d}x^{\mu}}{\mathrm{d}\zeta}\equiv p^{\mu}=\frac{\partial
H}{\partial\pi_{\mu}},\quad\frac{\mathrm{d}\pi_{\mu}}{\mathrm{d}\zeta}=-\frac{\partial
H}{\partial x^{\mu}}.$ (13)

The affine parameter $\zeta$ of the particle is related to its proper time
$\tau$ by the relation $\zeta=\tau/m$.

Due to the symmetries of the Schwarzschild spacetime (2) and the parabolic MF
(4), one can easily find the conserved quantities that are the energy and the
axial angular momentum of the particle and can be expressed as

$\displaystyle E$ $\displaystyle=$
$\displaystyle-\pi_{t}=mf(r)\frac{\mathrm{d}t}{\mathrm{d}\tau},$ (14)
$\displaystyle L$ $\displaystyle=$
$\displaystyle\pi_{\phi}=m\,r^{2}\sin^{2}\theta\,\frac{\mathrm{d}\phi}{\mathrm{d}\tau}+\frac{qB}{2}r^{w}(1-|\cos\theta|).$
(15)

The dynamical equations for the charged particle motion in the Cartesian
coordinates can be found by the transformation equations

$x=r\cos\phi\,\sin\theta,\quad y=r\sin\phi\,\sin\theta,\quad z=r\cos\theta.$
(16)

Introducing for convenience the specific parameters, energy ${\cal{E}}$, axial
angular momentum ${\cal{L}}$, and magnetic parameter ${\cal{B}}$, by the
relations

${\cal{E}}=\frac{E}{m},\quad{\cal{L}}=\frac{L}{m},\quad{\cal{B}}=\frac{qB}{2m},$
(17)

one can rewrite the Hamiltonian (11) in the form

$H=\frac{1}{2}f(r)p_{r}^{2}+\frac{1}{2r^{2}}p_{\theta}^{2}+\frac{1}{2}\frac{m^{2}}{f(r)}\left[V_{\rm
eff}(r,\theta)-{\cal{E}}^{2}\right],$ (18)

where $V_{\rm eff}(r,\theta;{\cal{L}},{\cal{B}})$ denotes the effective
potential given by the relation

$V_{\rm eff}(r,\theta)\equiv
f(r)\left[1+\left(\frac{{\cal{L}}}{r\sin{\theta}}-{\cal{B}}\,r^{w}\,\frac{1-|\cos\theta|}{r\sin{\theta}}\right)^{2}\right].$
(19)

The terms in the parentheses correspond to the central force potential given
by the specific angular momentum ${\cal{L}}$, and EM potential energy given by
the magnetic parameter ${\cal{B}}$. The Hamiltonian (18) can be divided into
dynamical part $H_{\rm D}$ (first two terms containing dynamical momenta
$p_{r},p_{\theta}$) and potential part $H_{\rm P}$ (last term).

The effective potential (19) shows clear symmetry
$({\cal{L}},{\cal{B}})\leftrightarrow(-{\cal{L}},-{\cal{B}})$ that allows to
distinguish the following two situations

  * $-$

minus configuration, here ${\cal{L}}>0,{\cal{B}}<0$ (equivalent to
${\cal{L}}<0,{\cal{B}}>0$) – MF and angular momentum parameters have opposite
signs and the LF is attracting the charged particle to the $z$-axis, towards
the BH.

  * +

plus configuration, here ${\cal{L}}>0,{\cal{B}}>0$ (equivalent to
${\cal{L}}<0,{\cal{B}}<0$) – MF and angular momentum parameters have the same
signs and the LF is repulsive, acting outward the BH.

In this article, without loss of generality, we use the positive angular
momentum of a particle ${\cal{L}}>0$, while the magnetic parameter ${\cal{B}}$
can be both positive or negative.

Figure 4:  For repulsive LF (${\cal{B}}>0$), we plot the radial positions of
the inISCOs for charged particles in equatorial plane, with different MF
configurations $w=0,0.1,0.75,1$. The dotted vertical line shows the position
of the BH horizon, while the dot-dashed vertical line indicates the position
of the inISCO for a neutral particle.  Figure 5: Energy in effective potential
minima, for off-equatorial (${\cal{B}}<0$) and in-equatorial (${\cal{B}}>0$)
plane. We set the inclination parameter of the MF to $w=3/4$. The black thick
curve is plotted for neutral particles. Each point of the dashed curve
represents the position of offISCO (or inISCO), while the numbers on solid
curves indicate the values of magnetic parameter ${\cal{B}}$.

###  2.1 Off- and in-equatorial circular orbits

The charged particle motion is confined by the energetic boundary given by

${\cal{E}}^{2}=V_{\rm eff}(r,\theta;{\cal{L}},{\cal{B}}).$ (20)

Let us properly investigate the features of the effective potential (19)
represented in Fig. 2 that enables us to demonstrate the general properties of
the charged particle dynamics, avoiding the necessity to solve the equations
of motion. Since the $V_{\rm eff}(r,\theta)$ function is symmetrical above
$\theta\in(0,\pi/2)$ and below $\theta\in(\pi/2,\pi)$ equatorial plane, we can
focus on the region $\theta\in(0,\pi/2)$ only, and we can substitute
$|\cos\theta|$ with $\cos\theta$ function here. The effective potential as a
function of coordinates $x$-$z$ is plotted in Fig. 2.

The effective potential $V_{\rm eff}(r,\theta)$ function is zero at BH horizon
and limits are

$\lim_{x\to\infty}V_{\rm eff}(x,z)=1,\quad\lim_{z\to\infty}V_{\rm
eff}(x,z)=1+\frac{{\cal{L}}^{2}}{x^{2}},$ (21)

for all values of magnetic parameters ${\cal{B}}$ and $w\in(0,1)$. Only for
value $w=1$, we have in the case of first limit $\lim_{x\to\infty}V_{\rm
eff}(x,z)=1+{\cal{B}}^{2}$.

The stationary points of the effective potential $V_{\rm eff}(r,\theta)$
function, where maxima or minima can exist, are given by the equations

$\partial_{r}V_{\rm
eff}(r,\theta;{\cal{L}},{\cal{B}})=0,\quad\partial_{\theta}V_{\rm
eff}(r,\theta;{\cal{L}},{\cal{B}})=0.$ (22)

These equations for the region above the equatorial plane $\theta\in(0,\pi/2)$
are

$\displaystyle({\cal{L}}+F)[({\cal{L}}$ $\displaystyle+$ $\displaystyle
F)+(r-2)(w-1)F+{\cal{L}}]\csc^{2}\theta$ (23) $\displaystyle+$ $\displaystyle
r^{2}=0,$ $\displaystyle({\cal{L}}$ $\displaystyle+$ $\displaystyle
F)({\cal{L}}\cos\theta-F)=0,$ (24)

where we define the supplementary function

$F(r,\theta)={\cal{B}}r^{w}(\cos\theta-1).$ (25)

The term containing $|\cos\theta|$ in $V_{\rm eff}(r,\theta)$ function (19)
will create the discontinuity in the first derivative of $V_{\rm
eff}(r,\theta)$ with respect to $\theta$ in an equatorial plane
$\theta=\pi/2$, and such points must be also added when we would like to
examine $V_{\rm eff}(r,\theta)$ function extrema.

Figure 6: Charged particles trajectories (blue) around BH with parabolic MF,
representing two astrophysically relevant cases: when MF influence is
negligible ($|{{\cal{B}}}|{\ll}1$, first row), and when LF is dominating over
gravity ($|{{\cal{B}}}|{\gg}1$, second row). The gray circle represents the
BH, gray curves MF lines, the dashed curve denotes the energy boundary, and
the black dot the initial position of the particle.

There exist two classes of effective potential extrema solutions: off-
equatorial plane extrema for attractive LF (${\cal{B}}<0$), and in-equatorial
plane extrema for repulsive LF (${\cal{B}}>0$).

For MF parameter ${\cal{B}}<0$, we have minima located at off-equatorial plane

$\cos(\theta)=\frac{{\cal{B}}r^{w}}{{\cal{B}}r^{w}-{\cal{L}}},$ (26)

while radial coordinate $r$ can be calculated using angular momentum condition
${\cal{L}}={\cal{L}}_{\rm off}(r)$, where

${\cal{L}}_{\rm off}(r)\equiv\frac{Y(r)+\sqrt{Y(r)^{2}+4(r-3)r^{2}}}{2(r-3)}.$
(27)

The new function $Y(r)$ is given by

$Y(r)={\cal{B}}\,r^{w}\left(w\left(2-r\right)+2\left(r-3\right)\right).$ (28)

Off-equatorial plane minima are limited by the minimal distance from the BH –
so-called innermost stable circular orbits (offISCOs) given by the equation

${\cal{B}}{\cal{L}}\left(r\left(w(1-w)+2\right)+2(w-3)w\right)r^{w-1}-{\cal{L}}^{2}+2r=0,$
(29)

where the particle angular momentum ${\cal{L}}$ is given by Eq. (27). The
radial positions of stable circular orbits as well as offISCO situated in the
off-equatorial plane for different values of magnetic parameter ${\cal{B}}<0$
and declination parameter $w$ are depicted in Fig. 3. Both offISCO and off-
equatorial stable circular orbits get closer to the $z$-axis with the increase
of the magnetic parameter ${\cal{B}}$. The offISCOs are depending also on MF
configuration, been close to the BH for pure magnetic monopole $w=0$ than in
the case of parabolic MF $w=3/4$.

Existence of the off-equatorial circular orbits in the background of
magnetized BHs with external asymptotically uniform or dipole MFs Kov-Stu-
Kar:2008:CLAQG: ; Kov-etal:2010:CLAQG: ; Kov-Kop-Kar-Koj:CQG-2013  inspired
the idea of investigation of the possible existence of charged off-equatorial
toroidal structures orbiting around magnetized BH. One should note, that BHs
can not have their own long-lived dipolar fields due to the no-hair theorem
Lyu-Mck:2011:PhRvD: ; Bransgrove-etal:2021:PhRvL: . Radiation belts filled
with energetic charged particles exists around many different astrophysical
objects (mostly with dipole magnetosphere), so also in the case of BH one can
search for observational consequences. As it can be seen in Fig. 3, the off-
equatorial plane minima are located closer to the $z$-axis (and within the jet
funnel) with increasing MF magnitude ${\cal{B}}$ forming off-equatorial
torodial structures in low density jet funnel filled with charged particle
orbiting around BH on the off-equatorial circular orbits.

For repulsive LF (${\cal{B}}>0$), we have minima located in equatorial plane
$\theta=\pi/2$, while radial coordinate can be calculated using angular
momentum condition ${\cal{L}}={\cal{L}}_{\rm in}(r)$, where

${\cal{L}}_{\rm
in}(r)\equiv\frac{Y(r)+\sqrt{{\cal{B}}^{2}w^{2}(r-2)^{2}r^{2w}+4(r-3)r^{2}}}{2(r-3)}.$
(30)

Effective potential minima in equatorial plane are limited by in-equatorial
innermost stable circular orbit (inISCO), given by the equation

$\displaystyle{\cal{B}}^{2}w^{2}\left(r-2\right)\left(wr^{2}-5wr+6w-r\right)r^{2w}$
$\displaystyle-{\cal
B}w\left(w\left(r^{2}-5r+6\right)-2r^{2}+11r-18\right)r^{w}\times$
$\displaystyle\times\sqrt{{\cal
B}^{2}w^{2}\left(r-2\right)^{2}r^{2w}+4\left(r-3\right)r^{2}}$
$\displaystyle+2\left(r-6\right)\left(r-3\right)r^{2}$ $\displaystyle=0.$ (31)

The radial position of inISCO for different magnitudes of the MF ${\cal{B}}$
and various values of parameter $w$, situated in an equatorial plane are
depicted in Fig. 4. The inISCO radius is constant $r_{\rm inISCO}=6$ only in
the case of split magnetic monopole magnetosphere $w=0$ Kol-Bar-
Jur:2019:RAGtime: . For other values of declination parameter $w$ of MF, the
radii of inISCOs lie in the interval $2\leq{r_{\rm inISCO}}{\leq}6$. For large
MF magnitude, the inISCO radii to the BH horizon $r_{\rm inISCO}=2$, however,
when the MF is negligible, we obtain the limiting case $r_{\rm inISCO}=6$.
Similar behaviour can be observed for uniform MF Kol-Stu-Tur:2015:CLAQG: ;
San-Pan-San:2022:PhRvD: .

Figure 7: Charged particle trajectories (blue) around BH with parabolic
magnetosphere on orbits around the off- and in-equatorial plane minima. Gray
disk represents the BH, gray curves the MF lines, the dashed curve the energy
boundary and the black dot the initial position of the particle.

The graphical behaviour of the particle energy ${\cal{E}}(r)$ as function of
radial coordinate $r$ for circular orbits located in both off- and in-
equatorial planes is shown in Fig. 5. Before the offISCO (or inISCO) position,
the energy of stable orbits decreases with the increase of radial coordinate,
but after the offISCO (or inISCO) position, the energy increases monotonically
as the radial distance increases. As it can be seen from Fig. 5, charged
particle on the off-equatorial plane circular orbit has lower binding energy
than particle on circular orbit in the equatorial plane. One can conclude that
off-equatorial circular orbits are less stable against perturbation than the
circular orbits located in equatorial plane.

Limiting cases, when the MF is negligible ($|{\cal{B}}|{\ll}1$) or when the MF
is dominating over gravity ($|{\cal{B}}|{\gg}1$) are more astrophysically
prevalent - examples of particle orbits for these two cases are plotted in Fig
6. In order to see clearly the interplay between magnetic LF and gravitation
BH attraction leading to strongly non-linear and chaotic dynamics, one should
use magnetic parameter $|{\cal{B}}|\sim{1}$. The effect of attractive
(${\cal{B}}<0$) and repulsive (${\cal{B}}>0$) LF is demonstrated on particle
trajectories around off- and in-equatorial minima in Fig 7. Effective
potential minima can naturally serve as traps for particles and one could
expect higher particle concentrations there, but the situation in/close to the
equatorial plane must be described by a more complex approach than our single
charged particle dynamics, for example by GRMHD of PIC simulations since
Keplerian accretion disk formed by conductive plasma is expected to orbit BH
in equatorial plane.

Figure 8: Influence of RR on charged particle dynamics in parabolic ($w=0.75$)
magnetosphere of BH. The initial conditions for particle dynamics are the same
as in Figs. 6 and 7. The radiation parameter $k$ has been chosen
unrealistically large in order to clearly demonstrate the effect of RR force.
Due to RR the particle spirals down to the BH. The damping of oscillations can
be observed in the first row. Time scales related to RR depend on MF parameter
${\cal{B}}$ and specific charge, they are relevant for the strong MF, while
negligible for weak one.

##  3 Radiating charged particle trajectories

Figure 9:  Charged particle chaotic dynamics (blue curve) in parabolic MF
(gray parabolic lines) around BH (gray circle) can be well demonstrated by the
PS plots (last column). In the first row, we see complete regular motion
represented by smooth the Kolmogorov–Arnold–Moser (KAM) torus section. In the
second row, we observe that this KAM torus is already starting to break and
chaos is emerging with this orbit. In the third row some special (unstable)
periodic orbit is presented, with only single point in PS. In the fourth row
one can observe completely chaotic orbit, with points fully covering its phase
space, limited only by particle’s energy (dashed black curve). We use MF
parameters ${\cal{B}}=2$, $w=0.75$ for all orbits.  Figure 10:  PSs for
chaotic particle orbit without (blue first row) and with (red second row) the
influence of RR. The same trajectories as in Fig. 9 have been used. The RR is
a non-linear damping force acting on a charged particle which can increase the
chaoticity of a trajectory, but can also lead to the particle collapse to the
BH.

In many astrophysically relevant scenarios, one cannot neglect the effects of
RR due to the synchrotron radiation of charges in the vicinity of BHs, which
are believed to be immersed in an external MF. Equation of motion for a
charged particle undergoing RR force in curved spacetime is non-trivial
Poisson:2004:LRR: . However, for elementary and subatomic particles, the
equations can be simplified to the form Tur-etal:2018:APJ:

$\displaystyle\frac{Du^{\mu}}{d\tau}=\frac{q}{m}F^{\mu}_{\,\,\,\nu}u^{\nu}$
(32)
$\displaystyle+\,\frac{q\,k}{m}\left(F^{\alpha}_{\,\,\,\beta;\mu}u^{\beta}u^{\mu}+\frac{q}{m}\left(F^{\alpha}_{\,\,\,\beta}F^{\beta}_{\,\,\,\mu}+F_{\mu\nu}F^{\nu}_{\,\,\,\sigma}u^{\sigma}u^{\alpha}\right)u^{\mu}\right),$

where semicolon denotes the covariant coordinate derivative. The Eq. (32),
which is a covariant form of the Landau-Lifshitz equation Lan-Lif:1975:CTF: ,
is a habitual second order differential equation which satisfies the principle
of inertia and does not contain runaway solutions. The exact equation of
motion expressed for our particular parabolic magnetosphere are rather complex
and hence we do not present exact formulas here, however we solve them
numerically. The effect of RR is characterized by the radiation parameter

$k=\frac{2}{3}\frac{q^{2}}{m}.$ (33)

Detailed discussion and derivation of (32) can be found in Tur-etal:2018:APJ:
. One should note, that in our approach the charged particle is radiating due
to acceleration given by the LF which has many orders of magnitude stronger
effect than the so-called geodesic synchrotron radiation Bre-etal:1973:PRD: ;
Tor-Dol:2022:PRD: ; in our approach such radiation on geodesic given by the
tail integral Tur-etal:2018:APJ:  can be neglected. Detailed analysis of the
particular case of the motion of charged particles around Schwarzschild and
Kerr BHs immersed into an external asymptotically uniform MF was presented in
Tur-etal:2018:APJ: ; Kol-Tur-Stu:2021:PRD: .

Trajectories of radiating charged particles around Schwarzschild BH immersed
in parabolic MF are depicted in Fig. 8. For both attractive (${\cal{B}}<0$)
and repulsive (${\cal{B}}>0$) LFs, the charged particle will lose energy and
angular momenta due to the RR force, which will force it to spiral down to the
BH. For radiating particle trajectories plotted in Fig. 8, we use the same
initial conditions as for non radiating ones presented in Figs. 6 and 7. For
repulsive (${\cal{B}}>0$) LF, we can observe oscillation damping due to
reaction leading to particle’s temporary stabilization on circular orbit
(“orbital parking”). This “parking” effect can be observed only in repulsive
${\cal{B}}>0$ case and similar effect for repulsive LF has been already
observed in the case of BH in uniform MF Tur-etal:2018:APJ: . For parabolic
MF, the situation is different than for uniform one, here the circular parking
orbit is not stable and particle will be slowly in-spiraling down to the BH.
The RR force acts as a damping force - the oscillatory character of the motion
of charged particles changes because of the loss of energy and angular
momentum and radiating charged particles will always be captured by the BH.
The off-equatorial plane minima are unstable for radiating particles and
charged particles can stay oscillating around off-equatorial plane minima only
for a limited amount of time, which could still be quite long since RR is
relevant only for relatively strong MFs.

It is well known that the neutral test particle motion in Schwarzschild
background is regular and fully integrable. On the other hand, the charged
particle dynamics around magnetized BHs is generally chaotic Kov-Kop-Kar-
Koj:CQG-2013 ; Pan-Kol-Stu:2019:EPJC: . The motion of charged particles is
regular (integrable) only for magnetic monopole field configuration $w=1$ due
to monopole’s spherical symmetry Kol-Bar-Jur:2019:RAGtime: . Strongly non-
linear and chaotic trajectories can be seen in the case of parabolic MF
$w=3/4$, see lower row in Fig. 6. Still, in parabolic MF not all orbits are
chaotic: the trajectories corresponding to minima of effective potential or
with small epicyclic oscillations around circular orbits are regular Kol-Stu-
Tur:2015:CLAQG: ; Tur-Stu-Kol:2016:PHYSR4: . Such almost circular orbits are
very important from the astrophysical point of view as they govern the thin
(Keplerian) accretion disks and the toroidal fluid configurations.

Charged particle chaotic behaviour can be analysed with the help of the
Poincaré section (PS), presented in Fig. 9. From the distribution of the
plotted points in the PS, one can distinguish whether the motion is chaotic or
not. For regular motion, the plotted points form a closed curve in the
2-dimensional $r-p_{r}$ plane, while in the case of the chaotic motion, PS
points are distributed apparently randomly in the allowed region. When the
energy of the particle is small, tori of small harmonic oscillation are not
broken, see first row of Fig. 9. Such tori get broken when the energy of the
particle is increased, see second row of Fig. 9, and the chaotic motion
appears. The sea of chaos spreads in the PS with the increase of the particle
energy, see fourth row of Fig. 9. The process of chaos increase with energy is
sometimes interrupted by regular orbits, see example of unstable periodic
orbit in third row of Fig. 9.

The comparison between PS for orbits of particles with and without the RR
force is plotted in Fig. 10. The RR force acts as a damping force leading the
radiating particle towards the BH attractor, hence trajectories at PS points
are moving radially downwards in time, as one can see the visible accumulation
of points eventually leading towards $r=2,p_{r}=-1$. The radiation losses are
highly non-linear during inspiral towards BH and particle may experience
multiple resonances Muk-Kop-Luk:2022:XXX: , which is demonstrated as
structures in radiating particle PS. Unstable periodic orbit, third column in
Fig. 10, will get perturbed by RR and become chaotic. Already chaotic, fourth
column in Fig. 10, will remain chaotic.

Figure 11:  Example of ionization mechanism for neutral particles in Keplerian
accretion disk. The LF and RR force acting on electrons and protons (ions) are
significantly different. In the presented case, the proton (blue) is
oscillating around effective potential minima in the equatorial plane and
influenced by LF only slightly (${\cal{B}}=0.01$). On the other hand, the
electron (red) is strongly influenced by the LF (${\cal{B}}=-10$) in the same
MF, orbiting MF lines, while moving up, down, and around the BH.

##  4 Ionized Keplerian disks

Single charged particle dynamics is a toy model for plasma description, but
even with this simple tool one can test combined gravitational and EM effects
of BH magnetosphere keeping control on numerical scheme, which could be
clouded in more complex GRMHD or PIC numerical simulations. To demonstrate the
influence of the MF on processes in BH vicinity, we consider an ionization of
Keplerian accretion disk consisting of electrically neutral test particles
following circular geodesics and orbiting a Schwarzschild BH immersed in a
parabolic MF, assuming the symmetry of Keplerian disk plane almost coinciding
with the equatorial plane of the Schwarzschild geometry.

Realistic ionization model of the originally neutral matter has been applied
e.g. for the description of the magnetic Penrose process Par-Wag-Dad:APJ:1986:
; Tur-etal:2020:APJ: , where the original neutral particle 1 splits into two
oppositely charged particles – 2nd and 3rd, with charges $q_{2}$ and $q_{3}$,
respectively. We thus assume conservation of the electric charge
$0=q_{2}+q_{3}$, and the canonical momentum

$\pi_{\alpha(1)}=\pi_{\alpha(2)}+\pi_{\alpha(3)},$ (34)

where $\pi_{\alpha(i)}$, $i=\\{1,2,3\\}$ denotes the canonical momentum of
1st, 2nd and 3rd particles, respectively. Due to the charge conservation, the
momentum conservation Eq. (34) takes the form

$p_{\alpha(1)}=p_{\alpha(2)}+q_{2}A_{\alpha}+p_{\alpha(3)}+q_{3}A_{\alpha}=p_{\alpha(2)}+p_{\alpha(3)}.$
(35)

In many realistic scenarios one of the newly created charged particles is much
more massive than the other one, $m_{2}/m_{3}\gg\leavevmode\nobreak\ 1$, such
as in the case of a neutral atom ionization. Example of ionization mechanism
for electrically neutral particle orbiting Schwarzschild BH immersed in
parabolic MF is shown in Fig. 11. The magnitude of BH magnetosphere is set to
$B$, but due to difference in masses of electrons and protons (ions), the LF
(${\cal{B}}$ parameter) and RR forces ($k$ parameter) are very different.
Therefore more massive charged particle (proton or ion) takes almost all the
initial momentum of the original neutral particle, and the dynamical influence
of the lighter charged particle (electron) can be neglected, so that

$p_{\alpha(1)}\approx p_{\alpha(2)}\gg p_{\alpha(3)}.$ (36)

In another realistic scenario, we can consider the Keplerian disk created by
plasma modelled as a quasi-neutral soup of charged particles: electrons and
ions, which are orbiting around the central object in circular orbits. If the
disk is dense enough, the mean free path of charged particles is much shorter
in comparison to the length of the orbit around the central object. This means
that the charged particles orbit the BH together as a collective neutral body
along circular geodesics as the influence of the MF is irrelevant to the
motion, but at the edges of the disk its density decreases substantially and
the charged particle motion starts to be significantly influenced by the MF as
the mean free path becomes comparable to the orbital length.

Figure 12:  Evolution of thin Keplerian accretion disks around Schwarzschild
BH immersed into external parabolic MF. Initially, the particles from
accretion disk follow circular geodesics with slight inclination from the
equatorial plane ($\theta_{0}=1.5$), while the MF is directed along the
$z$-axis. We set the declination parameter of the parabolic MF lines to
$w=3/4$, radiation parameter to $k=10^{-5}$, and the values for the MF
parameter to ${\cal{B}}=\\{\pm 10,\pm 1,\pm 0.1,\pm 0.01,0\\}$. The black
dashed closed curves are the ISCOs for neutral test particles, indicating the
inner edge of the Keplerian accretion disk. In the absence of the MF
(${\cal{B}}=0$ case, middle figure), all the orbits remain in their circular
shape and only the inclined razor-thin disk can be seen. If a slight EM
interaction is switched-on (${\cal{B}}=\pm 0.01$ cases, middle row), the
charged particles creating a disk that is initially almost perpendicular to
the MF lines start to oscillate with epicyclic frequencies around the circular
orbit in both radial and vertical directions; the accretion disk becomes
slightly thicker. Further, increasing the MF parameter (${\cal{B}}=\pm 0.1,\pm
1$ cases), the motion of charged particles becomes chaotic and the accretion
disk is either destroyed or modified into a thick toroidal structure. For
${\cal{B}}=-1$ case, all the particles are captured by the BH, and complete
destruction of the Keplerian disk can be observed. The LF dominates the
particle motion when the parameter of the MF is large enough ($|{\cal{B}}|\geq
10$ cases). The charged particles spiral down and up along the MF lines, and
move around the BH in the counter-clockwise (${\cal{B}}<0$) and clockwise
(${\cal{B}}>0$) directions.  Figure 13:  PSs for particles from ionized
Keplerian disk, see Fig. 12. Two upper plots correspond to trajectories
without RR: the left one (${\cal{B}}=0.1$) consists from regular trajectories
forming closed loops, while the right one (${\cal{B}}=2$) from chaotic
trajectories, where the points are distributed apparently randomly; the grey
dotted vertical lines show the inISCO positions. In the two lower plots, the
RR is included and we can see that some trajectories are attracted to the BH.
All radiating particle trajectories may eventually end up in BH, however for
small ${\cal{B}}$ and $k$, or if the orbit is far away from the BH, the
accretion takes much longer time.

In both scenarios neutral particles are ionized, while their mechanical
momenta are conserved. This simple ionization model (36) has been studied
previously in the field of rotating Kerr BH with uniform MF Stu-Kol:2016:EPJC:
; Kop-Kar:2018:APJ: ; Tur-etal:2021:PRD: . Chaoticity of the motion of the
ionized matter orbiting Schwarzschild BHs has been studied in detail in Pan-
Kol-Stu:2019:EPJC: .

Parabolic MF around non-rotating BH is given by $A_{\phi}$ component of the EM
four-potential, hence particle’s energy ($p_{t}$) is conserved during the
ionization process (35). Particles from neutral Keplerian disk are on bounded
orbits ($-p_{t}=E\leq 1$), hence they can not escape to infinity after
ionization since this would require ($E>0$), see (19). Charged particles are
able to escape along open parabolic MF lines to infinity if they obtain
initial kick in vertical direction – such a scenario has been described in the
uniform MF case in Zah-Fro-Sch:2013:PRD: ; Shi-Kim-Chi:2014:PRD: ; Bau-
Veg:2021:CQG: , where the parameters of the initial kick in vertical direction
were calculated. For Kerr BH, the situation is different since here BH
rotation generates “electric” component $A^{t}$ even in the case of a static
MF configuration. The presence of the $A^{t}$ component leads to the energy
redistribution in the formula (34), which can lead the particle to escape from
the bounded Keplerian circular orbit of BH with relatively high energy due to
the so-called magnetic Penrose process Stu-Kol:2016:EPJC: ; 2019Univ….5..125T
; Tur-etal:2020:APJ: ; Tur-Kol-Stu:Symmetry:2022: .

A test particle moving on the inclined circular orbit has initial position
$x^{\alpha}$ and four-velocity $u_{\alpha}$, which can be written in the form

$\displaystyle x^{\alpha}$ $\displaystyle=$
$\displaystyle(t,r,\theta,\phi)=(0,r_{0},\theta_{0},0),$ (37) $\displaystyle
u_{\alpha}$ $\displaystyle=$
$\displaystyle(u_{t},u_{r},u_{\theta},u_{\phi})=({\cal{E}},0,0,{\cal{L}}).$
(38)

Specific axial angular momentum ${\cal{L}}$ and specific energy ${\cal{E}}$
for neutral test particles on the inclined circular orbits from Keplerian
accretion disk are given by Wald:1984:book:

${\cal{L}}_{(I)}=\frac{r_{0}\sin\theta_{0}}{\sqrt{r_{0}-3}},\quad{\cal{E}}_{(I)}=\frac{r_{0}-2}{\sqrt{r_{0}^{2}-3r_{0}}}.$
(39)

Due to a simple ionization condition (35), using the definitions of the
specific energy and specific axial angular momentum (39), one can write for
the specific axial angular momentum ${\cal{L}}$ and the specific energy
${\cal{E}}$ of the ionized test particle

${\cal{L}}_{(II)}={\cal{L}}_{(I)}+A_{\phi}(r_{0},\theta_{0}),\quad{\cal{E}}_{(II)}={\cal{E}}_{(I)}.$
(40)

Only the specific angular momentum ${\cal{L}}$ is changed during the
ionization, while the particle specific energy ${\cal{E}}$ remains constant.
Charged particles from the ionized Keplerian disk can either be captured by
BH, oscillate around circular orbit, or move chaotically along parabolic MF
lines.

The influence of the MF parameter ${\cal{B}}$ on the fate of the Keplerian
disk orbiting a Schwarzschild BH in a position corresponding to the symmetry
plane of the parabolic MF is demonstrated in Fig. 12. When the disk remains
neutral or the MF is missing (${\cal{B}}=0$ case), all the orbits remain
circular – original inclined razor thin disk survives. For small values of the
magnetic parameter (${\cal{B}}=\pm 0.01$ cases), the charged particles will be
unsettled from original circular orbits and ionized disk will start to
oscillate in both radial and vertical directions; the ionized disk becomes
slightly thick due to the epicyclic motion of the charged matter. When the MF
is slightly increased (${\cal{B}}=\pm 0.1$ cases), complete destruction of the
inner part of the Keplerian disk occurs in the attractive ${\cal{B}}=-0.1$
case, as all the orbits end at the BH. In the repulsive ${\cal{B}}=0.1$ case,
the regular epicyclic motion has a tendency to chaoticity and the ionized disk
becomes thick. Further increasing the MF parameter (${\cal{B}}=\pm 1$ cases)
the charged particle motion enters a fully chaotic regime causing a
substantial thickening of the disk. For large values of the MF parameter
(${\cal{B}}=\pm 10$ cases), the LF becomes the leading force of the particle
motion. LF in both attractive and repulsive cases leads to similar results:
the charged particles go up and down along MF lines realizing epicyclic motion
around BH in clockwise (${\cal{B}}>0$) or counter-clockwise (${\cal{B}}<0$)
directions.

The influence of the RR on the ionized Keplerian accretion disks orbiting
Schwarzschild BH with parabolic MF can be observed for realistic values of the
radiation parameter $k$ only in a very long integration time; in Fig. 13 we
have used $k=0.02$ to make RR effect more visible. The RR force reduces the
chaotic behaviour of the particle motion. Neutral disk ionization and
subsequent radiation losses do not lead to effective filling of the off-
equatorial plane minima (${\cal{B}}<0$ case) with charged particles – off-
equatorial minima have higher energy than equatorial one, and off-equatorial
minima are not attractors for radiating charged particles. Different processes
similar to Earth’s MF trapping charged particles from the solar wind in the
Van Allen radiative belts can also be considered to form off-equatorial plane
structures. It is worth mentioning that that the physics of pair cascades
around BHs depend heavily on RR and Compton scattering Levinson-
Cerutti:2018:A&A: ; Chen-Yuan:2020:ApJ: ; Crinquand-etal:2020:PhRvL: .

Figure 14:  Off-equatorial radial profiles of the frequencies of small
harmonic oscillations $\nu_{r}$, $\nu_{\theta}$, $\nu_{\phi}$, and $\nu_{L}$
of charged particle orbiting BH of the mass of $M=10M_{\odot}$ in external
parabolic MF for different values of the parameters ${\cal{B}}$ and $w$.
Figure 15: In-equatorial radial profiles of the frequencies of small harmonic
oscillations $\nu_{r}$, $\nu_{\phi}$, and $\nu_{L}$ of charged particle
orbiting BH of the mass of $M=10M_{\odot}$ in external parabolic MF for
different values of the parameters ${\cal{B}}$ and $w$.

##  5 Particle oscillations and observed frequencies

Ionization of particles from neutral accretion disk in the equatorial plane
may lead to the disk oscillations for $|{\cal{B}}|<1$. If a test particle is
slightly displaced from the equilibrium position, corresponding to a stable
circular orbit, the particle starts to oscillate around the equilibrium
position realizing thus epicyclic motion governed by the linear harmonic
oscillations. In order to calculate the epicyclic frequencies of small
perturbations to such orbits, we employ a linear stability analysis to compute
the radial and vertical frequencies using a general method Kol-Stu-
Tur:2015:CLAQG: . Let us introduce a vector of canonical variables
$y^{K}=(x^{i},p_{i})$, so that we rewrite the equations of motion derived from
the Hamiltonian (18) in the form

$\dot{y}^{K}=f^{K}(y^{N}),$ (41)

where dot denotes the derivative with respect to the proper time $\tau$. We
are interested in the solutions to the Eq. (41) near equilibrium position
$y_{0}^{K}$ corresponding to a circular orbit, where $f^{K}(y_{0}^{N})=0$,
except $N=3$. Linearizing $y^{K}=y_{0}^{K}+\zeta^{K}$, where
$\zeta^{K}/y_{0}^{N}<<1$ denotes a small deviation vector, one gets a set of
linear differential equations with constant coefficients

$\dot{\zeta}^{K}=\left(\frac{\partial f^{K}}{\partial
y^{N}}\bigg{\rvert}_{y=y_{0}}\right)\zeta^{K}+O(\zeta^{2}).$ (42)

We decompose the solutions to this system into eigenvalues and eigenvectors of
the Jacobian matrix $\partial f^{K}/\partial y^{N}$. The eigenvalues $\lambda$
specify the rate at which trajectories with small differences in initial
conditions separate. Thus, we evaluate Eq. (42) using the Hamiltonian (18)
with the circular orbit values $y_{0}$ and find a pair of eigenvalues that we
interpret as radial and vertical angular frequencies given by

$\displaystyle\omega^{2}_{\mathrm{r}}$ $\displaystyle=$
$\displaystyle\frac{1}{2g_{rr}}\left(\frac{\partial^{2}H_{\rm p}}{\partial
r^{2}}+\frac{\partial^{2}H_{\rm p}}{\partial\theta^{2}}-\sqrt{\chi}\right),$
(43) $\displaystyle\omega^{2}_{\mathrm{\theta}}$ $\displaystyle=$
$\displaystyle\frac{1}{2g_{\theta\theta}}\left(\frac{\partial^{2}H_{\rm
p}}{\partial r^{2}}+\frac{\partial^{2}H_{\rm
p}}{\partial\theta^{2}}+\sqrt{\chi}\right),$ (44)

where the unknown function $\chi$ takes the form

$\chi=\left(\frac{\partial^{2}H_{\rm p}}{\partial
r^{2}}\right)^{2}+\left(\frac{\partial^{2}H_{\rm
p}}{\partial\theta^{2}}\right)^{2}+4\left(\frac{\partial^{2}H_{\rm
p}}{\partial r\partial\theta}\right)^{2}-2\frac{\partial^{2}H_{\rm
p}}{\partial r^{2}}\frac{\partial^{2}H_{\rm p}}{\partial\theta^{2}}.$ (45)

There exists the third fundamental angular frequency of the epicyclic particle
motion, namely the Keplerian (axial) frequency $\omega_{\mathrm{\phi}}$, given
by

$\omega_{\mathrm{\phi}}\equiv
u^{\phi}=g^{\phi\phi}\left({\cal{L}}-qA_{\phi}\right).$ (46)

The radial and vertical angular frequencies (43-44) measured by a local
observer for the case of in- and off- equatorial circular orbits have too
complicated character when expressed explicitly. For $w=2$ and in equatorial
plane, Eqs. (43) and (46) reduce to the uniform MF configuration Kol-Stu-
Tur:2015:CLAQG: , while ${\cal{B}}=0$ leads to the neutral particle case
Wald:1984:book: . In contrast to the neutral case, there exists angular
frequency so-called Larmor frequency for charged particles, which is
associated with the MF only, and given by a relation

$\omega_{\rm L}=\frac{q}{m}|\textbf{B}|,$ (47)

where $|\textbf{B}|$ can be found by Eq. (8). Due to the form of the
Hamiltonian $H_{\rm P}$ (18) in parabolic MF, the first and second derivatives
of $H_{\rm P}$ with respect to $\theta$ does not exist at the equatorial
plane, hence the construction of $\omega_{\mathrm{\theta}}$ is not possible at
the equatorial plane. However, one can still perturb circular orbit in
vertical directions around the equatorial plane, so that the perturbation will
lead to the non-harmonic oscillations. Equations for the harmonic and non-
harmonic oscillators are given, respectively, as follows

$\ddot{\theta}+\omega_{\mathrm{\theta}}^{2}\,\theta=0,\qquad\ddot{\theta}+C_{\theta}\,{\rm
sgn}(\theta)=0,$ (48)

where $C_{\theta}$ is a constant. Large vertical perturbation from equatorial
plane lead to chaotic dynamics, see the second row in Fig. 6.

Figure 16:  Power Spectral Density for three $r(\tau)$ (blue), $\theta(\tau)$
(red) and $\phi(\tau)$ (green) coordinates along the particle’s trajectory.
Fundamental frequencies calculated analytically for circular orbit and
presented in Figs. 14 and 15 are plotted as vertical dashed lines. The top row
of figures corresponds to trajectories without RR $k=0$, see, also Figs. 6 and
7. Cases $3$ and $4$ represent regular trajectories around off- and in-
equatorial plane minima with distinguishable main peaks; The case $2$
corresponds to the chaotic trajectory. Second row of figures is plotted for
radiating particles, where the higher frequency peaks (overtones) are
flattened due to radiation damping and the spectrum is more smooth (noisy) in
higher frequencies.

The locally measured angular frequencies $\omega_{\mathrm{r}}$,
$\omega_{\mathrm{\theta}}$, $\omega_{\mathrm{\phi}}$ and $\omega_{\rm L}$ are
connected to the angular frequencies measured by the static distant observers,
$\Omega_{\beta}$, by the gravitational redshift transformation

$\Omega_{\beta}=\frac{\mathrm{d}{\diamond}_{\beta}}{\mathrm{d}t}=\omega_{\beta}\,\frac{\mathrm{d}\tau}{\mathrm{d}t}=\omega_{\beta}\,\frac{f(r)}{{\cal{E}}(r)},$
(49)

where $(f(r)/{\cal{E}}(r))$ is the redshift coefficient, given by the function
$f(r)$ and the particle specific energy at the circular orbit ${\cal{E}}(r)$.
If the fundamental frequencies of the small harmonic oscillations related to
the distant observers, $\Omega_{\beta}$, are expressed in the physical units,
their dimensionless form has to be extended by the factor $c^{3}/GM$. Then the
frequencies of the charged particle radial and latitudinal harmonic
oscillations measured by a distant observer in physical units are given by

$\nu_{\beta}=\frac{1}{2\pi}\frac{c^{3}}{GM}\,\Omega_{\beta}.$ (50)

The behaviour of the fundamental frequencies
$\omega_{\mathrm{r}},\omega_{\mathrm{\theta}},\omega_{\mathrm{\phi}}$ and
their ratios can help us to distinguish different shapes of charged particle
epicyclic orbits in the vicinity of a stable circular orbit.

According to the Newtonian theory of gravitation, all the frequencies are
equal, $\omega_{\mathrm{r}}=\omega_{\mathrm{\theta}}=\omega_{\mathrm{\phi}}$,
giving ellipse as the only possible bounded trajectory of a test particle
around a gravitating spherically symmetric body. For uncharged particles
moving around a Schwarzschild BH, the relation
$\omega_{\mathrm{r}}<\omega_{\mathrm{\theta}}=\omega_{\mathrm{\phi}}$ holds,
and there exists a periapsis shift for bounded elliptic-like trajectories
implying the effect of relativistic precession that increases with decreasing
radius of the orbit as the strong gravity region is entered Ste-
Vie:1999:PHYSRL: . For charged particle frequencies around BH in presence of
MF, one gets in general
$\omega_{\mathrm{r}}\neq\omega_{\mathrm{\theta}}\neq\omega_{\mathrm{\phi}}$,
see Kol-Stu-Tur:2015:CLAQG: ; Tur-Stu-Kol:2016:PHYSR4: ; Kol-Tur-
Stu:2017:EPJC: .

The radial $\nu_{\mathrm{r}}$, vertical $\nu_{\mathrm{\theta}}$, orbital
$\nu_{\mathrm{\phi}}$ and Larmor $\nu_{\rm L}$ frequencies measured by a
distant observer in off-equatorial plane (${\cal{B}}<0$) for different MF
configurations as a function of radial coordinate $r$ are plotted in Fig. 14.
We see that orbital $\nu_{\mathrm{\phi}}$ and vertical $\nu_{\mathrm{\theta}}$
frequencies coincide for all considered MF configurations. When the attractive
LF is weak, the Larmor frequency is smaller than the orbital frequency,
however, it coincides with the orbital frequency as the magnetic parameter
${\cal{B}}$ is increased. For parabolic MF ($w=3/4$), we can see that the
radial frequency $\nu_{\mathrm{r}}$ is changing dramatically with the
parameter ${\cal{B}}$, while the vertical $\nu_{\mathrm{\theta}}$ and orbital
$\nu_{\mathrm{\phi}}$ frequencies remain almost the same. For strong
attractive LF (${\cal{B}}<-1$), the charged particle moves relatively slowly
around BH at the off-equatorial plane minima, while oscillating relatively
rapidly around MF lines. The radial frequency is too small for the case of
paraboloidal MF with $w=1$, so the radial profiles start to appear only far
away from the BH, which is in concord with the results presented for offISCO
in Fig. 3.

Fundamental frequencies for the case of in-equatorial plane (${\cal{B}}>0$)
are plotted in Fig. 15. When the repulsive LF is weak, the Larmor frequency
$\nu_{\rm L}$ is the smallest among all the frequencies, while the radial
$\nu_{\mathrm{r}}$ and orbital $\nu_{\mathrm{\phi}}$ frequencies coincide.
However, when the repulsive LF is strong (${\cal{B}}>1$), the radial
$\nu_{\mathrm{r}}$ and Larmor $\nu_{\rm L}$ frequencies have different
profiles in strong gravity regime close to the BH, but the frequencies
$\nu_{\mathrm{r}}$ and $\nu_{\rm L}$ will come closer as the radial distance
increases.

The effect of RR on Fourier spectra (Power Spectral Density) calculated from
charged particle trajectories is shown in Fig. 16. Particles moving close to
circular orbits have main peaks at fundamental frequencies given by Eqs. (43)
- (46). The higher harmonics (overtones) are also shown. The RR damping force
does not change main frequency peaks, but it lowers the higher harmonics
making the spectrum more flat and diluted in high frequency range (see cases
$3$ and $4$ in Fig. 16 in comparison with the cases $3R$ and $4R$,
respectively). The RR makes the motion more irregular (chaotic) on short time
scales. For both attractive (${\cal{B}}<0$ case 3) and repulsive
(${\cal{B}}>0$ case 4) LFs, the orbital motion ($\phi$, green) is considerably
suppressed in high frequency range. We interpret this effect as follows: since
the azimuthal motion (in $\phi$ direction) in quasi-circular orbits has the
highest curvature, hence larger acceleration with respect to other directions
of the motion, this component of the motion radiates away faster. Due to RR
the trajectory becomes more smooth on short time scales, since the high
curvature modes are dumped. Similar effect of RR damping can also be seen for
chaotic trajectory (case 2) due to the same reason. RR increases the
chaoticity in already chaotic case only slightly, although detectable using
the tools from non-linear analysis Pan-Ada-Mar:2022:EPJST: .

##  6 Astrophysical estimates

In order to relate our results to realistic astrophysical scenarios, in this
section we provide estimates of the most relevant parameters of the discussed
model in astrophysical situations. The MF strengths close to the BH surface
can be estimated from the Eddington value (MF energy density coinciding with
the density of accreting plasma) Rees:1984:ARAA:

$B_{\rm Edd}\sim 10^{4}\left(\frac{M}{10^{9}M_{\odot}}\right)^{-1/2}{\rm Gs}.$
(51)

According to Daly:2019:APJ: ; Pio-etal:2020:arXiv: , the characteristic values
of the MFs near the stellar mass and SMBHs are comparable with the Eddington
value. In Table 1 we give the values of $B_{\rm Edd}$, orbital timescale
$\tau_{c}$ of the particle at the radius $r=6M$ for two representative BH
masses: stellar mass BH ${10}\leavevmode\nobreak\ M_{\odot}$ and SMBH
${10^{9}}\leavevmode\nobreak\ M_{\odot}$.

| BH mass | orbital time | $B_{\rm Edd}$
---|---|---|---
stellar-mass BH |  $10^{1}M_{\odot}$ | $10^{-3}\leavevmode\nobreak\ {\rm{s}}$ | $10^{8}\leavevmode\nobreak\ {\rm{Gs}}$
SMBH |  $10^{9}M_{\odot}$ | $10^{4}\leavevmode\nobreak\ {\rm{s}}$ | $10^{4}\leavevmode\nobreak\ {\rm{Gs}}$
Table 1: Typical size, inISCO orbital timescale, and Eddington MF for stellar and SMBHs, representing characteristic values of many BH candidates.  MF $B$ | ${\cal B}_{\rm electron}$ | ${\cal B}_{\rm proton}$ |  ${\cal B}_{\rm Fe+}$ | ${\cal B}_{\rm dust}$
---|---|---|---|---
${10\leavevmode\nobreak\ {M}_{\odot}}$ |  |  |  |
$B=10^{8}$ Gs | $10^{10}$ | $10^{7}$ | $10^{5}$ | $10^{-2}$
$B=10^{4}$ Gs | $10^{6}$ | $10^{3}$ | $10^{1}$ | $10^{-6}$
$B=10^{0}$ Gs | $10^{2}$ | $10^{-1}$ | $10^{-3}$ | $10^{-10}$
$B=10^{-5}$ Gs | $10^{-3}$ | $10^{-6}$ | $10^{-8}$ | $10^{-15}$
${10^{9}\leavevmode\nobreak\ {M}_{\odot}}$ |  |  |  |
$B=10^{4}$ Gs | $10^{14}$ | $10^{11}$ | $10^{9}$ | $10^{2}$
$B=10^{0}$ Gs | $10^{10}$ | $10^{7}$ | $10^{5}$ | $10^{-2}$
$B=10^{-5}$ Gs | $10^{5}$ | $10^{2}$ | $10^{0}$ | $10^{-7}$
Sgr A* (${4\times 10^{6}\leavevmode\nobreak\ {M}_{\odot}}$) |  |  |
$B=10$ Gs | $10^{9}$ | $10^{6}$ | $10^{4}$ | $10^{-3}$
Table 2:  Magnetic parameter ${\cal B}$, Eq. (52), for various MF magnitudes
and different charged particles: electron, proton, ionized iron atom (one
electron lost), dust grain (one electron lost, $m\sim 10^{-18}$ g). Stellar
mass $M={10\leavevmode\nobreak\ {M}_{\odot}}$, SMBH
$M={10^{9}\leavevmode\nobreak\ {M}_{\odot}}$, and the Galactic centre BH (Sgr
A*) cases are considered.

The relative ratio between gravitational and magnetic LFs is represented by
the parameter ${\cal B}$. Restoring the world constants in Eq. (17), the
dimensionless parameter ${\cal B}$, widely used in the present paper, takes
the form

${\cal B}=\frac{|q|BGM}{2m_{e}c^{4}},$ (52)

reflecting the relative influence of the gravitational and MFs on the charged
particle motion. As an example, let us estimate ${\cal B}$ for an electron
around Galactic centre SMBH Sgr A*, which is currently one of the best known
SMBH candidates. The equipartition value of the MF strength at the event
horizon scale of Sgr A* is usually estimated to be around $10-100$ Gs Johnson-
etal:2015:Science: ; Eatough-etal:2013:Natur: . Mass of the BH is measured to
be around $4\times 10^{6}M_{\odot}$ (see, e.g. Eckart-etal:2017:FOP: ). Thus,
we get for electrons at the Galactic centre

${\cal B}_{\rm SgrA^{*}(\leavevmode\nobreak\
e^{-})}\approx\frac{|e|BGM}{2m_{e}c^{4}}\approx 1.8\times 10^{10}.$ (53)

For protons, the values of ${\cal B}_{\rm SgrA^{*}(\leavevmode\nobreak\
p^{+})}$ is of the order of $10^{7}$. In Table 2 we show the values of ${\cal
B}$ for electron, proton, iron atom (without one electron), and a dust
particle (without one electron) for various values of the MF strength and the
three BH types (stellar mass, supermassive, and Galactic centre BH). Another
robust way of estimating the MF in nearby the horizon region could be from
observed relativistic jet, which is believed to be launched by the Blandford-
Znajek mechanism Ripperda-etal:2022:ApJ: .

B [Gs] |  ${\cal{B}}$ |  $k{\cal{B}}$ |  $k{\cal{B}}^{2}$
---|---|---|---
$10^{15.8}$ | $\sim 10^{19}$ | $\sim 1$ | $\sim 10^{19}$
$10^{12}$ | $\sim 10^{15}$ | $\sim 10^{-4}$ | $\sim 10^{11}$
$10^{8}$ | $\sim 10^{11}$ | $\sim 10^{-8}$ | $\sim 10^{3}$
$10^{6.4}$ | $\sim 10^{9}$ | $\sim 10^{-9}$ | $\sim 1$
$10^{4}$ | $\sim 10^{7}$ | $\sim 10^{-12}$ | $\sim 10^{-5}$
$10^{0}$ | $\sim 10^{3}$ | $\sim 10^{-16}$ | $\sim 10^{-13}$
$10^{-2.9}$ | $\sim 1$ | $\sim 10^{-19}$ | $\sim 10^{-19}$
$10^{-5}$ | $\sim 10^{-2}$ | $\sim 10^{-21}$ | $\sim 10^{-23}$
Table 3: Magnitudes of different forces acting on radiating charged particle in the dimensionless form of Eq. (32) for different values of magnetic field strengths, compared to the gravitational force of the order of $\sim 1$. The estimates correspond to a relativistic electron in vicinity of a stellar mass BH with $M=10\leavevmode\nobreak\ {M}_{\odot}$.  B [Gs] |  $\tau_{\rm RR}^{e-}$ [s[ |  $\tau_{\rm RR}^{p+}$ [s] |  $\tau_{\rm RR}^{Fe+}$ [s]
---|---|---|---
$10^{12}$ | $10^{-16}$ | $10^{-6}$ | $10^{-1}$
$10^{8}$ | $10^{-8}$ | $10^{2}$ | $10^{7}$
$10^{4}$ | $1$ | $10^{10}$ | $10^{15}$
$10^{0}$ | $10^{8}$ | $10^{18}$ | $10^{23}$
$10^{-5}$ | $10^{18}$ | $10^{28}$ | $10^{33}$
Table 4: Typical cooling times of electron $\tau_{e}$, proton $\tau_{p}$, and
ionized iron atom (without one electron) for different MF strengths. For
charged dust particles, the synchrotron cooling is irrelevant.

RR force acting on a charged particle is represented by the dimensionless
parameter $k$ which has the form

$k=\frac{2}{3}\frac{q^{2}}{mGM}.$ (54)

The value of parameter $k$ is much lower than that of ${\cal B}$. For electron
orbiting stellar mass and SMBH we have, respectively

$\displaystyle k_{\rm BH}\sim 10^{-19}\quad{\rm for}\,\,\,M=10M_{\odot},$ (55)
$\displaystyle k_{\rm SMBH}\sim 10^{-27}\quad{\rm
for}\,\,\,M=10^{9}M_{\odot}.$ (56)

For electron around Sgr A*, one gets $k_{\rm SgrA^{*}}\sim 10^{-25}$. For
protons, the values of $k$ parameter is lower by the factor of
$m_{p}/m_{e}\approx 1836$, as in the case with ${\cal B}$. Despite the
weakness of parameter $k$ as compared to ${\cal B}$, it enters into the
equations for ultrarelativistic particles as $k{\cal B}^{2}$, which can make
the effect of the RR force considerably large. These values for RR have been
estimated for Sgr A* and M87 in Ripperda-etal:2022:ApJ:  and in more details
for M87 in Hakobyan-etal:2023:ApJ: .

The motion of radiating charged particle is described by Eq. (32), containing
several terms of different magnitudes. Let us estimate the contribution of
each of the term in the equation separately and compare the individual terms
with the gravitational force. The left hand side of (32) is purely
gravitational term given by the Christoffel symbols, which can be taken in the
dimensionless units to be $\sim 1$. On the right hand side of (32) one can
distinguish the LF $\sim{\cal{B}}$ and the contribution due to the RR
consisting from the terms, which are proportional to $\sim k{\cal{B}}$ and
$\sim k^{2}{\cal{B}}^{2}$. Representative values for all four terms forces are
given in Table 3. Here we assume the particles to be relativistic $v\sim 1$,
which allows us to use the magnitude estimates for distinguishing different
forces of particle dynamics, although the exact values of individual forces at
given point may differ. Large values of ${\cal B}$ in astrophysical settings
suggest that the effects of EM interaction and RR can not be neglected since
the LF and RR can be much stronger than gravity, but RR can overcome LF only
for particles with large specific charge (electron) and for very strong MF of
the order of $B>10^{15}$ Gs.

The RR timescale can be estimated from the following equation Tur-
etal:2018:APJ:

$\tau_{\rm RR}\approx\frac{3}{2}\frac{m^{3}c^{5}}{f(r)q^{4}B^{2}},$ (57)

where $f(r)$ is a lapse function. Typical synchrotron RR timescales of
particles in various MFs are given in Table 4. For the Galactic centre SMBH,
the electron cooling timescale is of the order of $10^{4}$ s, while orbiting
timescale at inISCO is $\sim 10^{3}$ s. In general, for ions the cooling
timescale is shorter than in case of electrons by the factor of
$(m_{p}/m_{e})^{3}\sim 10^{10}$. As one can see from the Table 4, the
synchrotron radiation losses of ions and heavier particles are relevant in the
presence of relatively strong MFs, i.e., much stronger than those
corresponding to BHs, e.g., in the vicinity of neutron stars.

| B [Gs] |  $\omega_{\rm L}$ [Hz]
---|---|---
inter stellar medium | $10^{-5}$ | $10^{2}$
Earth’s MF | $10^{0}$ | $10^{6}$
Galactic center SMBH | $10^{1}$ | $10^{8}$
typical AGN | $10^{4}$ | $10^{11}$
neutron star | $10^{12}$ | $10^{19}$
Table 5: Magnitudes of some relevant astrophysical MFs and their related
Larmor frequencies.

Fundamental frequencies for oscillating charged particle derived analytically
in previous section can be relevant in astrophysical settings, e.g. in
relation to the modulation of intensity flux (light curve) in quasi-periodic
oscillations Kol-Stu-Tur:2015:CLAQG: ; Kol-Tur-Stu:2017:EPJC: , observed in
many BH systems or in the model of the Galactic centre flares components Tur-
etal:2020:Flares: . Fundamental frequencies can also be closely related to
frequencies of radiated EM radiation Bre-etal:1973:PRD: , hence they manifest
themselves in produced EM spectrum and EM polarization. For strong MFs, the
spectrum will be dominated by the Larmor frequency, estimated in Table 5 for
some typical astrophysical MFs. As one can see that for typical AGNs ($10^{4}$
Gs) the Larmor frequency is of the order of hundreds of GHz, similar to the
observational frequency of the Event Horizon Telescope (230 GHz).

The synchrotron radiation losses of charged particles can be relevant in a
plasma surrounding BH, if the cooling timescale is shorter or comparable with
the mean free time between collisions of particles in plasma. The later can be
estimated from the equation 1981phki.book…..L

$\tau_{ee}\approx\frac{m_{e}^{2}(3k_{B}T_{e}/m_{e})^{3/2}}{\pi e^{4}n_{e}8{\rm
ln}\Lambda},$ (58)

where ${\rm ln}\Lambda$ is a Coulomb logarithm, given by

${\rm ln}\Lambda={\rm ln}\frac{3}{2e^{3}}\frac{(k_{B}T_{e})^{3/2}}{\sqrt{\pi
n_{e}}}.$ (59)

Note that $\tau_{ee}$ represents the collisional timescale between electron-
electron collisions in an electron-proton plasma. For proton-proton and
electron-proton collisions at the same temperature, the timescales are
typically longer. The representative estimates of collision timescales for
given densities and temperatures of plasma are given in Table 6.

It is interesting to estimate the electron collision timescales for the
environment of the Galactic centre SMBH. A particle density in a plasma
surrounding Sgr A* can be obtained from the observed flares properties Tur-
etal:2020:Flares: , estimated to be of the order of $n\sim 10^{7\pm 1}$cm-3.
Modeling the emission of simultaneous NIR and X-ray flares from Sgr A* the
Lorentz factor of electrons are estimated to be of the order of
$\gamma_{e}\sim 10^{3}$, which corresponds to the plasma with electron
temperature $T_{e}\sim 10^{9}$K. Then, the electron-electron collision
timescale is $\tau_{ee}$(Sgr A*)$=10^{4\pm 1}$s. Note that this estimate is
based on the flare observations, i.e., considered denser regions of the
Galactic centre compared to its average surrounding. Such estimate depends on
a number of assumptions - it may very well be that the density in a flaring
region is lower, as long as the magnetization is high enough to power
particles to Lorentz factors of $\sim 10^{3}$ Zhdankin-etal:2023:XXX . At
larger spatial scales, averaging the plasma environment, the collision
timescale is less by a few orders of magnitude 2018MNRAS.480.4408Z . Thus, the
plasma surrounding Galactic centre can be considered as collisionless, so the
considerations and estimates given in the current paper are fully applicable.
In principle, there is more magnetic energy available in dilute (high sigma)
regions, so flares could be powered from there. There are simulations showing
that flares from Sgr A* can be powered from magnetized regions, see Ripperda-
etal:2022:ApJ: ; Zhdankin-etal:2023:XXX . This however does not change the
main conclusion that the plasma is collisionless.

$n$ [cm-3] |  $T$ [K] |  $\tau_{ee}$ [s] |  $\tau_{ii}$ [s]
---|---|---|---
$10^{14}$ | $10^{10}$ | $10^{-1}$ | $1$
$10^{14}$ | $10^{6}$ | $10^{-7}$ | $10^{-5}$
$10^{10}$ | $10^{10}$ | $10^{3}$ | $10^{4}$
$10^{10}$ | $10^{6}$ | $10^{-3}$ | $10^{-1}$
$10^{6}$ | $10^{10}$ | $10^{7}$ | $10^{8}$
$10^{6}$ | $10^{6}$ | $10$ | $10^{3}$
Table 6: Collision timescale of electron-electron $\tau_{ee}$ and ion-ion
$\tau_{ee}$ collisions in a plasma of given particle density and temperature.
The plasma is considered quasi-neutral with equal temperatures of the
components. The values are rounded to the orders of magnitude.

##  7 Conclusions

The effects of MFs around BHs cannot be neglected in many astrophysical
settings. If the specific charge of a test particle $q/m$ is large enough, as
in the case of elementary particles, even a weak MF can significantly
influence the particle dynamics, the position of the inner edge of an ionized
Keplerian accretion disk, or the epicyclic frequencies for motion around BH.
In this article we have tested a simple but realistic model for BH
magnetosphere given by the parabolic MF and explored influence of RR damping
force on chaotic charge particle dynamics. We also briefly discussed the
astrophysical relevance of the presented theoretical model.

We studied the parabolic MF configuration around BH in general form for
different values of the inclination parameter $w$, limiting to the magnetic
split-monopole solution for $w=0$ as a special case. The parabolic
configuration of MF is supported by observations and GRMHD numerical
experiments Nak-etal:2018:APJ: . Generation of such MF geometry would require
the presence of electric current sheets at the equatorial plane provided by
ongoing accretion processes. The accretion disk is located in equatorial plane
and outside of the equatorial plane, the matter density is expected to be much
lower. Therefore, in discrete particle approach used in the present paper, the
existence of charged particle off-equatorial orbits become especially
important.

We have analysed stability of circular orbits of charged particles and found
that the off-equatorial stable circular orbits may exist only in the case of
attractive LF (${\cal{B}}<0$), while the stable circular orbits at the
equatorial plane can be observed only in the case of repulsive LF
(${\cal{B}}>0$). The effective potential minima situated at the off-equatorial
plane get closer to the $z$-axis when the MF becomes stronger. In the case of
the in-equatorial plane, the radius of inISCO of the charged particle
decreases as the MF strength increases, approaching the BH horizon $r_{\rm
inISCO}=2$ at large values of the magnetic parameter ${\cal{B}}$. In the case
of the split magnetic monopole, the inISCO (or offISCO) remains at a constant
radius $r_{\rm inISCO}=6$.

We examined the regular and chaotic behaviour of the motion of charged
particles with the help of PSs, which is a good tool illustrating the
difference between regular and chaotic dynamics. We have found that charged
particle motion around BH in parabolic MF can be chaotic under certain initial
conditions. This result is in contrast with the neutral particle case or the
case with the monopole MF configuration, where due to the symmetry of the EM
and gravitational fields, chaos does not occur and the motion is always
regular Kol-Bar-Jur:2019:RAGtime: . We have also compared the change of the
chaoticity of trajectories when RR force in included in the equations of
motion and found that the RR acts as a damping force, creating new attractors
in phase space, which are situated on circular orbit at the equatorial plane
of the BH. We have found that the off-equatorial plane orbits are always
unstable when the radiation is taken into account, since the RR force tends to
shift the charged particle orbits to their guiding centers, which coincides
with the position of the BH. The RR force has damping character of the
influence on the charged particle dynamics and the final stage of the
radiating particle motion is eventually a collapse into the BH if the
integration time is long enough. The particle loses its energy and angular
momentum due to RR cooling, which leads to spiraling down to the BH. For
positive value of magnetic parameter ${\cal{B}}>0$, an oscillating charged
particle experiencing the RR force radiates its oscillations away and settles
down to the quasi-circular orbits, before eventually falling down into the BH.

We have also studied the influence of the RR force on the structure of a thin
Keplerian accretion disk orbiting Schwarzschild BH in external parabolic MF.
The effect of RR can be observed even for a relatively small values of
${\cal{B}}$. For $|{\cal{B}}|\gg 1$, one can observe oscillations in the disk
dynamics, for $|{\cal{B}}|\sim 1$, the disk trajectories enter into the
turbulent chaotic regime and for $|{\cal{B}}|\ll 1$, the charged particle
trajectories stick to MF lines. Due to RR force, disk loses energy and angular
momentum and particle will end up into BH. The time scale of the destruction
of the Keplerian disk by radiation strongly depends on the strength of MF and
radiation parameter $k{\cal{B}}$.

We have calculated the fundamental frequencies of a charged particle and
produced EM spectrum, which can be useful in determining time scales of
various astrophysical processes around BH and fitting observed data. Test
particle frequencies are fundamental in the sense that they will appear in any
more advanced and complicated model for BH neighborhood description, like
GRMHD simulation as a zero approximation Mis-Klu-Fra:2019:MNRAS: . It is
worthwhile mentioning that our analytical solutions are very important to test
the codes that evolve particles (e.g. GRPIC), for example in the codes of
Parfrey-etal:2019:PhRvL: ; Bacc-etal:2019:ApJS: ; Bacchini-etal:2020:ApJS: .
One example where particle frequencies could be crucial is the quasi-periodic
oscillations of the X-ray power density observed in quasars and microquasars.
Fitting these quasi-periodic oscillation frequencies in timing spectra can
help to verify BH magnetosphere model Kol-Tur-Stu:2017:EPJC: ; Sha-
etal:2021:EPJC: . Another observational window for BH accretion processes
could be non-thermal synchrotron emission produced by charged particles moving
in the BH magnetosphere Fro-etal:2021:ARXIV: . One charged particle dynamics
and its oscillations around effective potential minima could be used to
construct synthetic charged particle synchrotron spectra Sho:2015:PRD: ;
Crinquand-etal:2022:PhRvL: ; Hakobyan-etal:2023:ApJ:

We have also estimated the most relevant parameters of the discussed model in
astrophysical situations and found the limits on the relevance of the effects
of the LF and RR for different particle types and BH systems. For that, we
gave the estimates of orbital, radiation, and collisional timescales of
particles. We have shown that applied to the Galactic centre BH the plasma
environment can be considered as collisionless and test particle approach can
be used there. We also discussed the role of MF in fundamental frequencies,
observable in relation to the so-called quasi-periodic oscillations.

Among the possible extensions of the presented parabolic magnetosphere model,
which could be interesting to explore in future are the calculation of the
radiated EM spectra and polarization, and inclusion of the BH rotation into
the model. In relation to the existence of the off-equatorial orbits found in
this paper, it would also be interesting to explore, whether the stable off-
equatorial structures, existing in discrete particle approach can be found in
numerical GRMHD or PIC simulations. This is a prediction that needs to be
tested.

## Acknowledgments

The work is supported by the Research Centre for Theoretical Physics and
Astrophysics, Institute of Physics, Silesian University in Opava and the GAČR
23-07043S project.

## References

  * (1) R.P. Eatough, H. Falcke, R. Karuppusamy, K.J. Lee, D.J. Champion, E.F. Keane, G. Desvignes, D.H.F.M. Schnitzeler, L.G. Spitler, M. Kramer, B. Klein, C. Bassa, G.C. Bower, A. Brunthaler, I. Cognard, A.T. Deller, P.B. Demorest, P.C.C. Freire, A. Kraus, A.G. Lyne, A. Noutsos, B. Stappers, N. Wex, Nature501, 391 (2013).  DOI 10.1038/nature12499
  * (2) R.A. Daly, The Astrophysical Journal886(1), 37 (2019).  DOI 10.3847/1538-4357/ab35e6
  * (3) R. Gold, J.C. McKinney, M.D. Johnson, S.S. Doeleman, The Astrophysical Journal837(2), 180 (2017).  DOI 10.3847/1538-4357/aa6193
  * (4) Event Horizon Telescope Collaboration, K. Akiyama, J.C. Algaba, A. Alberdi, W. Alef, R. Anantua, K. Asada, R. Azulay, A.K. Baczko, D. Ball, M. Baloković, J. Barrett, B.A. Benson, D. Bintley, L. Blackburn, R. Blundell, W. Boland, K.L. Bouman, G.C. Bower, H. Boyce, M. Bremer, C.D. Brinkerink, R. Brissenden, S. Britzen, A.E. Broderick, D. Broguiere, T. Bronzwaer, D.Y. Byun, J.E. Carlstrom, A. Chael, C.k. Chan, S. Chatterjee, K. Chatterjee, M.T. Chen, Y. Chen, P.M. Chesler, I. Cho, P. Christian, J.E. Conway, J.M. Cordes, T.M. Crawford, G.B. Crew, A. Cruz-Osorio, Y. Cui, J. Davelaar, M. De Laurentis, R. Deane, J. Dempsey, G. Desvignes, J. Dexter, S.S. Doeleman, R.P. Eatough, H. Falcke, J. Farah, V.L. Fish, E. Fomalont, H.A. Ford, R. Fraga-Encinas, P. Friberg, C.M. Fromm, A. Fuentes, P. Galison, C.F. Gammie, R. García, Z. Gelles, O. Gentaz, B. Georgiev, C. Goddi, R. Gold, J.L. Gómez, A.I. Gómez-Ruiz, M. Gu, M. Gurwell, K. Hada, D. Haggard, M.H. Hecht, R. Hesper, E. Himwich, L.C. Ho, P. Ho, M. Honma, C.W.L. Huang, L. Huang, D.H. Hughes, S. Ikeda, M. Inoue, S. Issaoun, D.J. James, B.T. Jannuzi, M. Janssen, B. Jeter, W. Jiang, A. Jimenez-Rosales, M.D. Johnson, S. Jorstad, T. Jung, M. Karami, R. Karuppusamy, T. Kawashima, G.K. Keating, M. Kettenis, D.J. Kim, J.Y. Kim, J. Kim, J. Kim, M. Kino, J.Y. Koay, Y. Kofuji, P.M. Koch, S. Koyama, M. Kramer, C. Kramer, T.P. Krichbaum, C.Y. Kuo, T.R. Lauer, S.S. Lee, A. Levis, Y.R. Li, Z. Li, M. Lindqvist, R. Lico, G. Lindahl, J. Liu, K. Liu, E. Liuzzo, W.P. Lo, A.P. Lobanov, L. Loinard, C. Lonsdale, R.S. Lu, N.R. MacDonald, J. Mao, N. Marchili, S. Markoff, D.P. Marrone, A.P. Marscher, I. Martí-Vidal, S. Matsushita, L.D. Matthews, L. Medeiros, K.M. Menten, I. Mizuno, Y. Mizuno, J.M. Moran, K. Moriyama, M. Moscibrodzka, C. Müller, G. Musoke, A. Mus Mejías, D. Michalik, A. Nadolski, H. Nagai, N.M. Nagar, M. Nakamura, R. Narayan, G. Narayanan, I. Natarajan, A. Nathanail, J. Neilsen, R. Neri, C. Ni, A. Noutsos, M.A. Nowak, H. Okino, H. Olivares, G.N. Ortiz-León, T. Oyama, F. Özel, D.C.M. Palumbo, J. Park, N. Patel, U.L. Pen, D.W. Pesce, V. Piétu, R. Plambeck, A. PopStefanija, O. Porth, F.M. Pötzl, B. Prather, J.A. Preciado-López, D. Psaltis, H.Y. Pu, V. Ramakrishnan, R. Rao, M.G. Rawlings, A.W. Raymond, L. Rezzolla, A. Ricarte, B. Ripperda, F. Roelofs, A. Rogers, E. Ros, M. Rose, A. Roshanineshat, H. Rottmann, A.L. Roy, C. Ruszczyk, K.L.J. Rygl, S. Sánchez, D. Sánchez-Arguelles, M. Sasada, T. Savolainen, F.P. Schloerb, K.F. Schuster, L. Shao, Z. Shen, D. Small, B.W. Sohn, J. SooHoo, H. Sun, F. Tazaki, A.J. Tetarenko, P. Tiede, R.P.J. Tilanus, M. Titus, K. Toma, P. Torne, T. Trent, E. Traianou, S. Trippe, I. van Bemmel, H.J. van Langevelde, D.R. van Rossum, J. Wagner, D. Ward-Thompson, J. Wardle, J. Weintroub, N. Wex, R. Wharton, M. Wielgus, G.N. Wong, Q. Wu, D. Yoon, A. Young, K. Young, Z. Younsi, F. Yuan, Y.F. Yuan, J.A. Zensus, G.Y. Zhao, S.S. Zhao, Astrophysical Journal Letters910(1), L13 (2021).  DOI 10.3847/2041-8213/abe4de
  * (5) Z. Stuchlík, M. Kološ, J. Kovář, P. Slaný, A. Tursunov, Universe 6(2), 26 (2020).  DOI 10.3390/universe6020026
  * (6) R.D. Blandford, R.L. Znajek, Monthly Notices of the Royal Astronomical Society179, 433 (1977)
  * (7) D.L. Meier, _Black Hole Astrophysics: The Engine Paradigm_ (Springer Berlin, Heidelberg, 2012)
  * (8) A. Galishnikova, A. Philippov, E. Quataert, F. Bacchini, K. Parfrey, B. Ripperda, arXiv e-prints arXiv:2212.02583 (2022).  DOI 10.48550/arXiv.2212.02583
  * (9) M. Nakamura, K. Asada, K. Hada, H. Pu, S. Noble, C. Tseng, K. Toma, M. Kino, H. Nagai, K. Takahashi, J.C. Algaba, M. Orienti, K. Akiyama, A. Doi, G. Giovannini, M. Giroletti, M. Honma, S. Koyama, R. Lico, K. Niinuma, F. Tazaki, The Astrophysical Journal868(2), 146 (2018).  DOI 10.3847/1538-4357/aaeb2d
  * (10) O. Porth, Event Horizon Telescope Collaboration, The Astrophysical Journal Supplement243(2), 26 (2019).  DOI 10.3847/1538-4365/ab29fd
  * (11) M. Kološ, A. Janiuk, in _RAGtime 20-22: Workshops on black holes and neutron stars_ , ed. by S. Hledík, Z. Stuchlík (2020), pp. 153–164
  * (12) B. Crinquand, B. Cerutti, G. Dubus, K. Parfrey, A. Philippov, Astronomy and Astrophysics650, A163 (2021).  DOI 10.1051/0004-6361/202040158
  * (13) J.C. McKinney, R. Narayan, Monthly Notices of the Royal Astronomical Society375(2), 513 (2007).  DOI 10.1111/j.1365-2966.2006.11301.x
  * (14) J.C. McKinney, R. Narayan, Monthly Notices of the Royal Astronomical Society375(2), 531 (2007).  DOI 10.1111/j.1365-2966.2006.11220.x
  * (15) A. Tchekhovskoy, R. Narayan, J.C. McKinney, The Astrophysical Journal711(1), 50 (2010).  DOI 10.1088/0004-637X/711/1/50
  * (16) B. Turimov, Z. Stuchlík, J. Rayimbayev, A. Abdujabbarov, Phys. Rev. D103(12), 124039 (2021).  DOI 10.1103/PhysRevD.103.124039
  * (17) B. Ripperda, M. Liska, K. Chatterjee, G. Musoke, A.A. Philippov, S.B. Markoff, A. Tchekhovskoy, Z. Younsi, Astrophysical Journal Letters924(2), L32 (2022).  DOI 10.3847/2041-8213/ac46a1
  * (18) K. Hirotani, R. Krasnopolsky, H. Shang, K.i. Nishikawa, M. Watson, The Astrophysical Journal908(1), 88 (2021).  DOI 10.3847/1538-4357/abd3a6
  * (19) J. Kovář, Z. Stuchlík, V. Karas, Classical and Quantum Gravity 25(9), 095011 (2008).  DOI 10.1088/0264-9381/25/9/095011
  * (20) J. Kovář, O. Kopáček, V. Karas, Z. Stuchlík, Classical and Quantum Gravity 27(13), 135006 (2010).  DOI 10.1088/0264-9381/27/13/135006
  * (21) J. Kovář, O. Kopáček, V. Karas, Y. Kojima, Classical and Quantum Gravity 30, 025010 (2013)
  * (22) M. Lyutikov, J.C. McKinney, Phys. Rev. D84(8), 084019 (2011).  DOI 10.1103/PhysRevD.84.084019
  * (23) A. Bransgrove, B. Ripperda, A. Philippov, Physical Review Letters 127(5), 055101 (2021).  DOI 10.1103/PhysRevLett.127.055101
  * (24) M. Kološ, D. Bardiev, B. Juraev, in _Proceedings of RAGtime 20-21: Workshops on black holes and neutron stars_ (2019), pp. 145–152
  * (25) M. Kološ, Z. Stuchlík, A. Tursunov, Classical and Quantum Gravity 32(16), 165009 (2015).  DOI 10.1088/0264-9381/32/16/165009
  * (26) S. Shaymatov, P. Sheoran, S. Siwach, Phys. Rev. D105(10), 104059 (2022).  DOI 10.1103/PhysRevD.105.104059
  * (27) E. Poisson, Living Reviews in Relativity 7, 6 (2004).  DOI 10.12942/lrr-2004-6
  * (28) A. Tursunov, M. Kološ, Z. Stuchlík, D.V. Galtsov, The Astrophysical Journal861, 2 (2018).  DOI 10.3847/1538-4357/aac7c5
  * (29) L.D. Landau, E.M. Lifshitz, _The Classical Theory of Fields_ (Oxford: Pergamon Press, 1976)
  * (30) R.A. Breuer, P.L. Chrzanowksi, H.G. Hughes, C.W. Misner, Phys. Rev. D 8, 4309 (1973).  DOI 10.1103/PhysRevD.8.4309.  URL https://link.aps.org/doi/10.1103/PhysRevD.8.4309
  * (31) T. Torres, S.R. Dolan, Phys. Rev. D 106, 024024 (2022).  DOI 10.1103/PhysRevD.106.024024.  URL https://link.aps.org/doi/10.1103/PhysRevD.106.024024
  * (32) M. Kološ, A. Tursunov, Z. Stuchlík, Phys. Rev. D103(2), 024021 (2021).  DOI 10.1103/PhysRevD.103.024021
  * (33) R. Pánis, M. Kološ, Z. Stuchlík, arXiv e-prints arXiv:1905.01186 (2019)
  * (34) A. Tursunov, Z. Stuchlík, M. Kološ, Phys. Rev. D93(8), 084012 (2016).  DOI http://dx.doi.org/10.1103/PhysRevD.93.084012
  * (35) S. Mukherjee, O. Kopacek, G. Lukes-Gerakopoulos, arXiv e-prints arXiv:2206.10302 (2022)
  * (36) S. Parthasarathy, S.M. Wagh, S.V. Dhurandhar, N. Dadhich, The Astrophysical Journal307, 38 (1986).  DOI 10.1086/164390
  * (37) A. Tursunov, Z. Stuchlík, M. Kološ, N. Dadhich, B. Ahmedov, The Astrophysical Journal895(1), 14 (2020).  DOI 10.3847/1538-4357/ab8ae9
  * (38) Z. Stuchlík, M. Kološ, European Physical Journal C 76, 32 (2016).  DOI 10.1140/epjc/s10052-015-3862-2
  * (39) O. Kopáček, V. Karas, The Astrophysical Journal853, 53 (2018).  DOI 10.3847/1538-4357/aaa45f
  * (40) A. Tursunov, B. Juraev, Z. Stuchlík, M. Kološ, Phys. Rev. D104(8), 084099 (2021).  DOI 10.1103/PhysRevD.104.084099
  * (41) A.M. Al Zahrani, V.P. Frolov, A.A. Shoom, Phys. Rev. D87(8), 084043 (2013).  DOI 10.1103/PhysRevD.87.084043
  * (42) R. Shiose, M. Kimura, T. Chiba, Phys. Rev. D90(12), 124016 (2014).  DOI 10.1103/PhysRevD.90.124016
  * (43) J. Bautista, I. Vega, Classical and Quantum Gravity 38(15), 155016 (2021).  DOI 10.1088/1361-6382/ac0e19
  * (44) A. Tursunov, N. Dadhich, Universe 5(5), 125 (2019).  DOI 10.3390/universe5050125
  * (45) A. Tursunov, M. Kološ, Z. Stuchlík, Symmetry 14(3), 482 (2022).  DOI 10.3390/sym14030482
  * (46) R.M. Wald, _General relativity_ (University of Chicago Press, Chicago, 1984)
  * (47) A. Levinson, B. Cerutti, Astronomy and Astrophysics616, A184 (2018).  DOI 10.1051/0004-6361/201832915
  * (48) A.Y. Chen, Y. Yuan, The Astrophysical Journal895(2), 121 (2020).  DOI 10.3847/1538-4357/ab8c46
  * (49) B. Crinquand, B. Cerutti, A. Philippov, K. Parfrey, G. Dubus, Physical Review Letters 124(14), 145101 (2020).  DOI 10.1103/PhysRevLett.124.145101
  * (50) L. Stella, M. Vietri, Physical Review Letters 82, 17 (1999).  DOI 10.1103/PhysRevLett.82.17
  * (51) M. Kološ, A. Tursunov, Z. Stuchlík, European Physical Journal C 77, 860 (2017).  DOI 10.1140/epjc/s10052-017-5431-3
  * (52) R. Pánis, K. Adámek, N. Marwan, European Physical Journal Special Topics (2022).  DOI 10.1140/epjs/s11734-022-00686-4
  * (53) M.J. Rees, Annual Review of Astronomy and Astrophysics22, 471 (1984).  DOI 10.1146/annurev.aa.22.090184.002351
  * (54) M.Y. Piotrovich, A.G. Mikhailov, S.D. Buliga, T.M. Natsvlishvili, arXiv e-prints arXiv:2004.07075 (2020)
  * (55) M.D. Johnson, V.L. Fish, S.S. Doeleman, D.P. Marrone, R.L. Plambeck, J.F.C. Wardle, K. Akiyama, K. Asada, C. Beaudoin, L. Blackburn, R. Blundell, G.C. Bower, C. Brinkerink, A.E. Broderick, R. Cappallo, A.A. Chael, G.B. Crew, J. Dexter, M. Dexter, R. Freund, P. Friberg, R. Gold, M.A. Gurwell, P.T.P. Ho, M. Honma, M. Inoue, M. Kosowsky, T.P. Krichbaum, J. Lamb, A. Loeb, R.S. Lu, D. MacMahon, J.C. McKinney, J.M. Moran, R. Narayan, R.A. Primiani, D. Psaltis, A.E.E. Rogers, K. Rosenfeld, J. SooHoo, R.P.J. Tilanus, M. Titus, L. Vertatschitsch, J. Weintroub, M. Wright, K.H. Young, J.A. Zensus, L.M. Ziurys, Science 350, 1242 (2015).  DOI 10.1126/science.aac7087
  * (56) A. Eckart, A. Hüttemann, C. Kiefer, S. Britzen, M. Zajaček, C. Lämmerzahl, M. Stöckler, M. Valencia-S, V. Karas, M. García-Marín, Foundations of Physics 47, 553 (2017).  DOI 10.1007/s10701-017-0079-2
  * (57) H. Hakobyan, B. Ripperda, A.A. Philippov, Astrophysical Journal Letters943(2), L29 (2023).  DOI 10.3847/2041-8213/acb264
  * (58) A. Tursunov, M. Zajaček, A. Eckart, M. Kološ, S. Britzen, Z. Stuchlík, B. Czerny, V. Karas, The Astrophysical Journal897(1), 99 (2020).  DOI 10.3847/1538-4357/ab980e
  * (59) E.M. Lifshitz, L.P. Pitaevskii, _Physical kinetics_ (Pergamon Press, 1981)
  * (60) V. Zhdankin, B. Ripperda, A.A. Philippov, arXiv e-prints arXiv:2302.05276 (2023).  DOI 10.48550/arXiv.2302.05276
  * (61) M. Zajaček, A. Tursunov, A. Eckart, S. Britzen, Monthly Notices of the Royal Astronomical Society480(4), 4408 (2018).  DOI 10.1093/mnras/sty2182
  * (62) B. Mishra, W. Kluźniak, P.C. Fragile, Monthly Notices of the Royal Astronomical Society483(4), 4811 (2019).  DOI 10.1093/mnras/sty3124
  * (63) K. Parfrey, A. Philippov, B. Cerutti, Physical Review Letters 122(3), 035101 (2019).  DOI 10.1103/PhysRevLett.122.035101
  * (64) F. Bacchini, B. Ripperda, O. Porth, L. Sironi, The Astrophysical Journal Supplement240(2), 40 (2019).  DOI 10.3847/1538-4365/aafcb3
  * (65) F. Bacchini, B. Ripperda, A.A. Philippov, K. Parfrey, The Astrophysical Journal Supplement251(1), 10 (2020).  DOI 10.3847/1538-4365/abb604
  * (66) M. Shahzadi, M. Kološ, Z. Stuchlík, Y. Habib, European Physical Journal C 81(12), 1067 (2021).  DOI 10.1140/epjc/s10052-021-09868-1
  * (67) C.M. Fromm, A. Cruz-Osorio, Y. Mizuno, A. Nathanail, Z. Younsi, O. Porth, H. Olivares, J. Davelaar, H. Falcke, M. Kramer, L. Rezzolla, arXiv e-prints arXiv:2111.02518 (2021)
  * (68) A.A. Shoom, Phys. Rev. D92(12), 124066 (2015).  DOI 10.1103/PhysRevD.92.124066
  * (69) B. Crinquand, B. Cerutti, G. Dubus, K. Parfrey, A. Philippov, Physical Review Letters 129(20), 205101 (2022).  DOI 10.1103/PhysRevLett.129.205101