1. Introduction
Symplectic integrators [
1,
2,
3,
4,
5] are geometric integrators that ensure the long-term stability of energy errors. By constructing a symmetric symplectic algorithm, truncation errors of odd orders can be effectively eliminated. Symplectic algorithms can be categorized into explicit symplectic algorithms and implicit symplectic algorithms. Explicit symplectic integrators have been widely used in black hole spacetime. The constructions of explicit symplectic integrators in Schwarzschild black holes and Reissner–Nordström black holes are very simple because their Hamiltonians can be easily split into four and five independently integrable parts, respectively [
6,
7]. Then, Wang et al. [
8] applied this method to a Hamiltonian for the description of charged particles moving around a Reissner–Nordström–(anti)–de Sitter black hole with an external magnetic field and successfully constructed second- and fourth-order explicit symplectic integrators by splitting that Hamiltonian into six analytically solvable pieces. They discovered that these integrators exhibit exceptional long-term behavior in preserving the boundedness of Hamiltonian errors, irrespective of whether the orbits are ordered or chaotic, as long as appropriate step sizes are employed. For Kerr spacetime, in order to utilize an explicit symplectic integrator, a time transformation function needs to be applied to the Hamiltonian of Kerr geometry [
9]. On this basis, Huang et al. [
10] successfully constructed time-transformed explicit symplectic schemes in a magnetized deformed Schwarzschild black hole spacetime, and these schemes had good performance in stabilizing energy errors without secular drift. Recently, Zhang et al. [
11] successfully constructed a time-transformed explicit symplectic integrator for a Hamiltonian system that describes the movement of charged particles around a non-Schwarzschild black hole that is surrounded by an external magnetic field.
As a subject of continuous exploration [
12,
13,
14], the phenomenon of chaos is gradually becoming known. Until now, although there has been a clear understanding of the chaos characteristics of Newtonian systems, the chaos characteristics of relativistic gravity systems are quite vague. For example, a double fixed-center problem is integrable in Newtonian mechanics, but it has strong chaotic behavior in the relativistic framework [
15,
16,
17,
18]. Therefore, the study of chaotic dynamics has been a hot topic in recent years [
19,
20,
21,
22,
23,
24,
25,
26,
27,
28]. In the study conducted by Gonçalo A. S. Dias and José P. S. Lemos [
29], they successfully constructed thin-shell electrically charged wormholes in
d-dimensional general relativity with a cosmological constant. The energy conditions of the static wormhole solutions were analyzed by using the Darmois–Israel formalism, and a linearized stability analysis was carried out. According to the study of A. M. Al Zahrani et al. [
30], they discovered that the trajectories of charged particles moving around a weakly magnetized static nonrotating black hole have three different types of asymptotic behavior, and the escape of a particle to infinity exhibits features similar to those of a diffusion process. In the research of Cao et al. [
31], it was found that the presence of a fourth constant of motion that is similar to the Carter constant leads to regular behavior of charged particles in moving around a Kerr–Newman black hole that is surrounded by cloud strings, quintessence, and an electromagnetic field. Stuchlík, Z. et al. [
32] presented a comprehensive review on the impact of cosmic repulsion and external magnetic fields on the accretion disks surrounding rotating black holes, as well as the jets associated with these rotating configurations. Their study revealed that in ionized Keplerian disks, an external magnetic field combined with the chaotic scattering process crucially influences the fate of the disks, enabling protons and ions with energy reaching
eV to exist simultaneously in the vicinity of rotating supermassive black holes with
while immersed in a sufficiently strong magnetic field with
.
As the study of chaos progresses, numerous methods have been developed to identify chaos. The most typical ways are the Poincaré section method, the Lyapunov index [
33], the local Lyapunov index and its spectral distribution [
34], the fast Lyapunove index [
35,
36], the smaller alignment index [
37], and the generalized alignment index [
38]. Each method has its own unique strengths and weaknesses. For the Poincaré section, its advantage is that it directly reflects the dynamic properties of the system, while its limitation is that it is only suited for studying conservative systems with two degrees of freedom. However, the Lyapunov index measures the intensity of chaos by studying how the orbital deviation vectors change over time, and it can reflect the average separation ratio of two adjacent orbitals over time. This method is applicable to systems of any dimension, but it often requires a long calculation time. On the other hand, the fast Lyapunov index, as the name suggests, can identify chaos faster than the traditional Lyapunov index can. The smaller alignment index, which can accurately identify the possible alignment of two multidimensional vectors, is a fast and effective method for detecting chaos. The generalized alignment index describes the linear independence of deviation vectors, and it is a natural extension of the smaller alignment index for more than two deviation vectors. However, its calculation is considerably complex.
In this study, we construct an efficient explicit symplectic algorithm in a non-rotating Konoplya and Zhidenko black hole and employ it to calculate some chaos indicators. We found the most effective approach and applied it to examine the dynamics of charged particles moving around the black hole. The rest of this study is structured as follows. In
Section 2, we construct an explicit symplectic scheme (
) of the Hamiltonian in the non-rotating Konoplya and Zhidenko black hole spacetime and an explicit scheme of the Hamilton variational equations. In
Section 3, the effectiveness of the SALI is verified, and the influences of different algorithms on the effectiveness of the SALI are investigated. In
Section 4, we examine the influence of the parameter
.
Section 5 is devoted to conclusions and discussions.
3. Numerical Evaluations
We consider the parameters
,
,
Q = 8.9
,
, and the proper time step
. The initial conditions are
,
, and
. The initial value of
(positive) is determined by Equation (
9). As shown in
Figure 1, we obtain the Hamiltonian energy errors of different numerical algorithms. We can see that the Hamiltonian error of
is between
and
, while that of
is between
and
. Furthermore, these errors are bounded and stable. In contrast, the error of
increases with the calculation time. This is a natural result because
is not a geometric integrator, which means that it does not preserve the symplectic structure of Hamiltonian systems.
Then, we utilize some methods to determine the orbital dynamics of the charged test particles moving around the black hole. They are the Poincaré section, FLI, and SALI.
The definition of the FLI is expressed as [
43]
where
and
represent the Euclidean norms of the deviation vector
at times
and
, respectively. If the orbit is chaotic, the FLI will grow exponentially. Conversely, if it is regular, the FLI will exhibit linear growth.
The SALI is given by [
37]
where
and
, are unit vectors. When the orbit is chaotic, the two unit vectors tend to align with each other, leading to SALI
. However, for a regular orbit, the values of the SALI will exhibit bounded fluctuations around a positive number.
Figure 2a illustrates the time evolution of the FLIs for four different orbits, while
Figure 2b depicts the corresponding Poincaré sections. We can see in
Figure 2a that when
, the FLI increases exponentially over time, indicating that the corresponding orbit exhibits strong chaos; for the weaker chaotic orbit
, the FLI grows more slowly. For ordered orbits
and
, the FLIs present linear growth. The Poincaré sections, which intuitively indicate whether an orbit is chaotic, are shown in
Figure 2b. In the case of chaotic orbits, the points in the Poincaré sections are randomly distributed throughout the plane. In contrast, for ordered orbits, the Poincaré sections display one or more closed curves.
In
Figure 3, we use different algorithms to obtain the SALIs of the system, whose initial conditions and parameters are the same as those in
Figure 2. In
Figure 3a, we adopt
to calculate the equations of motion and
to calculate the variational equations (
–
). The results indicate that the SALI of the strongly chaotic orbit with
approaches zero, while the SALI of the weakly chaotic orbit with
appears to be more stable than the former. For ordered orbits with
and
, the SALIs fluctuate around two constants. A calculation time of
is needed here, and this is faster than the FLI method in
Figure 2a.
In
Figure 3b, we utilize the
–
method to calculate the SALIs, and the results are nearly identical to those in
Figure 3a. In
Figure 3c,d, the SALIs are obtained by using global symplectic algorithms
–
and the
–three-particle method (due to the SALI requiring the computation of two deviation vectors and both deviation vectors needing triple integration by simultaneously using
, this is referred to as the three-particle method), respectively. For the global symplectic method, only a calculation time of
is needed to distinguish between chaos and order, making it the fastest approach. However, due to the rapid variation in the deviation vectors resulting from this calculation method, it becomes difficult to differentiate weakly chaotic orbits from strongly chaotic orbits. On the other hand, the SALIs obtained from the
–three-particle method lose their effectiveness.
Next, we are interested in examining how the parameters
E,
Q,
L, and
influence the orbital dynamics. Due to the computational advantage of SALIs in terms of the required calculation time, we employ the SALI and the Poincaré section to investigate the impacts of the four parameters on the dynamic transition from order to chaos, and the results are shown in
Figure 4. In
Figure 4(a1), it can be observed that if
or
, the SALIs are greater than
, which is the boundary between chaos and order in [
10], indicating that the orbits are ordered. However, for
and
, the corresponding orbits are chaotic. In
Figure 4(a2), when
and
, the points in the Poincaré sections form two rings, which are known as the Komogorov–Arnold–Moser (KAM) torus. This represents that these two orbits are orderly. When
and
, the orbits appear to be clearly chaotic.
By changing the value of parameter
, we find that, in
Figure 4(b1), the orbits are regular when
and
. The orbits transition into chaos when
. In addition, if
exceeds 9.1, the SALI is smaller than
, and as the value of
increases, the SALI decreases, and the extent of chaos is strengthened. In
Figure 4(b2), the motions of test particles are chaotic in
and
and regular in
and
.
We now investigate whether the angular momentum of the test particles influences the dynamic transition from order to chaos. In
Figure 4(c1), strong chaos is observed when
. When
, although there are some weakly chaotic orbits, most of the trajectories exhibit ordered behavior. In
Figure 4(c2), the values of angular momentum corresponding to the four orbits are, respectively,
,
,
, and
. Notably, the orbit with
is a chaotic orbit, while the orbit with
is ordered. If we set
, the orbit seems to be a typical weak chaotic orbit, and for
, the points in the figure form a closed ring, which indicates an ordered orbit.
Finally, we explore the transition of orbits when increasing the value of the magnetic field parameter
Q. In
Figure 4(d1), it is observed that the orbits transition from ordered to chaotic when
and back to ordered when
. The Poincaré section corresponding to
forms a closed torus, indicating a regular orbit. On the other hand, the Poincaré section of
becomes a complicated KAM torus, which is characterized by the presence of six small loops. These small loops belong to the same trajectory and form a chain of islands [
13]. Such a torus is considered regular, but it can easily lead to the onset of resonance and chaos. The orbit corresponding to
is a weakly chaotic orbit, while the orbit for
is strongly chaotic. However, the strength of the chaos in the latter is weaker than that depicted in
Figure 4(a2).
In order to conduct a more thorough investigation of the mixing of chaotic and regular motion, we focus on the exploration of chaotic particle dynamics in the
X–
Z plane.
Figure 5 illustrates the chaotic nature of particle trajectories in certain
X–
Y planes. The step size for
r is 0.05, while the step size for
is
. Deeper shades of red are utilized to illustrate regions in which particles exhibit more prominent chaotic behavior, while lighter shades of blue are employed to indicate areas with less pronounced chaotic behavior. We can clearly observe that when the value of
r falls within the range of approximately
to
, the particles demonstrate a significant manifestation of chaotic behavior. However, when
R exceeds
, the behavior of the particles becomes orderly and organized. It is noticeable that in certain areas, the chaotic behavior of particles sharing the same value of
r changes as
varies. The most prominent occurrence of this mixing of the chaotic and regular motion can be observed in the transitional region between chaos and order (roughly, when r falls within the range of 8.7 to 9.8).
4. The Effect of
In this section, we examine the impact of on effective potentials and stable circular orbits and its influence on chaos dynamics through a comparison with the magnetized Schwarzschild black hole.
From Equations (
8) and (
9), we have
where
V represents the effective potential on the equatorial plane
, which characterizes the radial movement of a charged particle within an accretion disk encircling the black hole:
The relationship between the effective potential and
can be observed in
Figure 6, which shows that as
increases, the effective potential decreases—see, in particular, the black line with
. For the circular motion of the particles, the following conditions must be satisfied:
The presence of local extrema in the effective potential indicates the existence of circular orbits. If the second derivative of the effective potential is
, the potential has local minimal values. These values correspond to the existence of stable circular orbits (SCOs). Specifically, when
, such a stable circular orbit is referred to as the innermost stable circular orbit (ISCO). The radii of some SCOs and ISCOs in
Figure 6 are listed in
Table 1. To enhance the visual representation, we have drawn a rough sketch of the SCOs and ISCOs in
Figure 7 for three different values of
: 0, 4, and 8. By referring to
Table 1 and
Figure 7, it is evident that as the parameter
increases, the radius of the ISCO increases, while the radius of the SCO decreases.
By setting the value of
to zero and maintaining the same initial conditions as those depicted in
Figure 4, we obtain the correlations between the SALIs and the energy
E, charge
Q, and angular momentum
L in a magnetized Schwarzschild black hole, and these are displayed in
Figure 8. Comparing
Figure 8a with
Figure 4(a1), we can observe that the relationship between
E and the SALI in
Figure 8a is generally similar to that in (a1) in
Figure 4. By comparing
Figure 8b,c with (c1) and (d1) in
Figure 4, respectively, we can see that when
is set to zero, there is no chaotic behavior observed when the angular momentum
L and charge
Q are changed.
5. Conclusions and Discussion
In summary, compared with the traditional method, the fourth-order explicit symplectic scheme is a superior method in non-rotating Konoplya and Zhidenko black hole spacetime. It ensures that the energy errors remain below and do not exhibit growth over time. The numerical results demonstrate that the SALI calculated with the global symplectic algorithm is the fastest method for distinguishing chaotic orbits. However, it is challenging to distinguish weak chaos when using this method. On the other hand, the – and – methods are capable of distinguishing weak chaos, but a longer calculation time is required.
The SALI is employed to study how the parameters affect the dynamic transition from order to chaos. The overall results are dynamic properties of all orbital changes as these parameters change. Specifically, when and the energy E go beyond a certain value, the orbits tend to remain chaotic. Conversely, when the angular momentum L exceeds a certain threshold, while there may still be some chaotic orbits, the majority of the orbits are regular. If we set the magnetic field parameter Q within a suitable range, a chaotic orbit can be easily found. However, if Q is below or exceeds a certain value, the orbits tend to be mainly regular.
By altering the initial conditions of the particles, a distribution plot illustrating the smaller Lyapunov exponent on the X–Z plane of the black hole was generated. It was observed that the particle orbits manifest an exceptionally diverse and intricate structure, and when the orbital radius of the particles r falls within the range of to , we can easily find the mixing of chaotic and regular motion.
The influence of the newly introduced parameter can cause a contraction in the SCO and an expansion in the ISCO. Furthermore, it can lead to chaotic behavior of the particles, which is not present in Schwarzschild black holes.
Researching the motion of charged particles around a black hole contributes to understanding the mechanisms of black hole accretion and provides insights into the initial formation process of accretion disks. Building upon this foundation, future considerations may involve studying the dynamic behavior of photons around a black hole, which can offer perspectives for investigating black hole shadows.