Higher-order geodesic deviations and orbital precession in a Kerr-Newman space-time

TL;DR

This paper introduces a higher-order geodesic deviation method to analyze orbital precession in Kerr-Newman space-times, allowing precise calculations for planetary orbits without relying on Newtonian approximations.

Contribution

The paper generalizes orbital precession calculations using higher-order geodesic deviations in Kerr-Newman space-times, avoiding Newtonian and post-Newtonian approximations for small eccentricities.

Findings

Accurate orbital trajectories without Newtonian approximations.

Generalization of perihelion precession in Kerr-Newman space-times.

Applicable for arbitrary values of GM/RC².

Abstract

A novel approximation method in studying the perihelion precession and planetary orbits in general relativity is to use geodesic deviation equations of first and high-orders, proposed by Kerner et.al. Using higher-order geodesic deviation approach, we generalize the calculation of orbital precession and the elliptical trajectory of neutral test particles to Kerr$-$Newman space-times. One of the advantage of this method is that, for small eccentricities, one obtains trajectories of planets without using Newtonian and post-Newtonian approximations for arbitrary values of quantity ${G M}/{R c^2}$.

Equations

$$\Delta \phi =\frac{6\pi GM}{a(1-e^2)},$$

$$\Delta \phi =\frac{6\pi GM}{a}(1+e^2+e^4+e^6+\cdots ).$$

$$\frac{D^2{n}^{\mu }}{D{\tau }^2}=R_{\lambda \nu \kappa }^{\mu }{u}^{\lambda }{n}^{\nu }{u}^{\kappa }.$$

$$\frac{d^2{n}^{\mu }}{d{\tau }^2}+2{\Gamma }_{\kappa \nu }^{\mu }{u}^{\kappa }\frac{d{n}^{\nu }}{d\tau }+{\partial }_{\nu }{\Gamma }_{\kappa \sigma }^{\mu }{u}^{\kappa }{u}^{\sigma }{n}^{\nu }=0.$$

$$x^{\mu }(\tau ,\lambda )=x^{\mu }(\tau ,{\lambda }_0)+(\lambda -{\lambda }_0)\frac{\partial x^{\mu }}{\partial \lambda }{|}_{(\tau ,{\lambda }_0)}+\frac{1}{2!}(\lambda -{\lambda }_0{)}^2\frac{{\partial }^2{x}^{\mu }}{\partial {\lambda }^2}{|}_{(\tau ,{\lambda }_0)}+\cdots .$$

$$b^{\mu }=\frac{D{n}^{\mu }}{D\lambda }=\frac{\partial n^{\mu }}{\partial \lambda }+{\Gamma }_{\lambda \nu }^{\mu }{n}^{\lambda }{n}^{\nu }.$$

$$x^{\mu }(\tau ,\lambda )=x^{\mu }(\tau ,{\lambda }_0)+ϵn^{\mu }+\frac{1}{2!}{ϵ}^2(b^{\mu }-{\Gamma }_{\lambda \nu }^{\mu }{n}^{\lambda }{n}^{\nu })+\cdots ,$$

$$x^{\mu }(\tau ,\lambda )=x_0^{\mu }(\tau )+\delta x^{\mu }(\tau )+\frac{1}{2}{\delta }^2{x}^{\mu }(\tau )+\cdots ,$$

$$\frac{D^2{b}^{\mu }}{D{s}^2}+R_{\rho \lambda \sigma }^{\mu }{u}^{\lambda }{b}^{\rho }{u}^{\sigma }=(R_{\rho \lambda \sigma ;\nu }^{\mu }-R_{\sigma \nu \rho ;\lambda }^{\mu })u^{\lambda }{u}^{\sigma }{n}^{\rho }{n}^{\nu }+4R_{\rho \lambda \sigma }^{\mu }{u}^{\lambda }\frac{D{n}^{\rho }}{Ds}{n}^{\sigma }.$$

$$\frac{d^2{b}^{\mu }}{d{\tau }^2}$$

$$d{s}^2$$

$${\rho }^2=r^2+a^2{\cos }^2\theta ,\Delta =r^2+a^2+Q^2-2Mr,$$

$$u^t=\dot{t}=\epsilon -\frac{(Q^2-2Mr)[(r^2+a^2)\epsilon -al]}{\Delta {\rho }^2},$$

$$u^{\phi }=\dot{\phi }=\frac{l}{{\rho }^2{\sin }^2\theta }-\frac{a[(Q^2-2Mr)\epsilon +al]}{\Delta {\rho }^2},$$

$${\rho }^2{\dot{r}}^2=[(r^2+a^2)\epsilon -al{]}^2-\Delta [r^2+(l-a\epsilon {)}^2],$$

$${\rho }^2\ddot{r}=2r\epsilon [(r^2+a^2)\epsilon -al]-r\Delta -(r-M)[r^2+(l-a\epsilon {)}^2].$$

$$u^r$$

$$u^t$$

$$u^{\phi }$$

$$\left(\begin{array}[]{cccc}\frac{d^{2}}{d\tau^{2}}&\alpha\frac{d}{d\tau}&0&0\\ \beta\frac{d}{d\tau}&\frac{d^{2}}{d\tau^{2}}+\kappa&0&\gamma\frac{d}{d\tau}\\ 0&0&\frac{d^{2}}{d\tau^{2}}+\delta&0\\ 0&\eta\frac{d}{d\tau}&0&\frac{d^{2}}{d\tau^{2}}\end{array}\right)\left(\begin{array}[]{cc}n^{t}\\ n^{r}\\ n^{\theta}\\ n^{\varphi}\end{array}\right)=\left(\begin{array}[]{cc}0\\ 0\\ 0\\ 0\end{array}\right),$$

$$\alpha =\frac{2(MR-Q^2)f_5}{R^3{f}_2\sqrt{f_1}},\beta =\frac{2(MR-Q^2)f_2}{R^3\sqrt{f_1}},\kappa =-\frac{3(MR-\frac{4}{3}{Q}^2)f_2}{R^4{f}_1},$$

$$\gamma =-\frac{2\sqrt{MR-Q^2}}{R}(1+\frac{a\sqrt{MR-Q^2}}{R^2})\frac{f2}{\sqrt{f1}},\delta =\frac{(MR-Q^2)f3}{R^{4}f1},\eta =\frac{2\sqrt{MR-Q^2}f4}{R^{3}f2\sqrt{f1}},$$

$$f_1$$

$$f_3$$

$$f_5$$

$${\omega }_{\theta }=\frac{\sqrt{\frac{M}{R}-\frac{Q^2}{R^2}}}{R}\sqrt{\frac{f_3}{f_1}},$$

$$n^{\theta }=-n_0^{\theta }\cos ({\omega }_{\theta }\tau ).$$

$$\frac{d^2{n}^r}{d{s}^2}+{\omega }^2{n}^r=0,$$

$$\omega =\frac{\sqrt{\frac{M}{R}-\frac{Q^2}{R^2}}}{R}\sqrt{\frac{f_6}{f_1}},$$

$$f_6=1-\frac{6M}{R}+\frac{8a}{R}\sqrt{\frac{M}{R}-\frac{Q^2}{R^2}}-\frac{3a^2}{R^2}+\frac{4Q^2}{R^2}+\frac{Q^2}{R^2}\left(1-\frac{M}{R}+\frac{a^2}{R^2}\right){\left(\frac{M}{R}-\frac{Q^2}{R^2}\right)}^{-1}.$$

$$n^r=-n_0^r\cos (\omega \tau ),$$

$$n^{\phi }$$

$$n^t$$

Full text

PACS numbers: 04.20.-q, 04.70. Bw, 04.80. Cc

Keywords: geodesic deviation, Kerr$-$Newman metric, perihelion

1 Introduction

Over the past century, the theory of General Relativity (GR) has paved the way for many predictions including, planet’s perihelion advance [1], bending of light when passing close by the sun, gravitational red shift [2], accelerating expansion of the universe [3] and the existence of gravitational waves [4], which have been confirmed by observational evidence [5].

The precession of planet’s perihelion is a famous problem in general relativity, see Figure 1. The residual perihelion advance of Mercury is 43 arc-seconds per century. Newtonian mechanics accounts for 531 arc-seconds per century by computing the perturbing influence of other planets in the solar system, while the predicted value of precession by observation is 574 arc-second per century. This residual value was obtained from observation by Urbain Le Verrier in 1859. To explain this discrepancy Li Verrier suggested a new planet (named Vulcon) near the sun which turned out not to exist [6].

Einstein [1] was the first to investigate the motion of planets as test particles moving along geodesics in a Schwarzschild space-time which, as is well known, leads to elliptic integrals. By expanding the integrand in powers of the small parameter $\frac{GM}{r}$ and keeping only linear terms in the expansion, he found the well-known formula for perihelion advance [7]

[TABLE]

where $G$ is the gravitational constant, $M$ is the mass of the central object, $a$ and $e$ are the major-semi axis and the eccentricity of the ellipse respectively. For small eccentricities the above relation can be expanded as a power series

[TABLE]

Interestingly, even for planet Mercury with greatest orbital eccentricity, the series converges rapidly such that the inclusion of the second term, that is $(1+e^{2})$, will lead to a result that differs only by 0.18% from the result predicted by relation (1), which is below the actual error bars. This would suggest that if orbital trajectories of test particles can be presented in a power series of eccentricity, rapid conversion and precise predictions will be achieved. This is what will be done in what follows, using the geodesic deviation method.

A novel approximation method for calculating the perihelion precession and finding trajectory of test particles in various gravitational fields without introducing the Newtonian and Post-Newtonian approximation has been suggested by Kerner et.al. [8]. In two previous papers, the authors calculate the perihelion precession and trajectory of test particles in Schwarzschild and Kerr space-times [9]–[13]. Our aim is to generalize the geodesic deviation method to the Kerr-Newmann geometry. The motion of neutral and charged test particles along closed orbits in Kerr-Newman metric has been analysed in an exhaustive manner in many papers [14]–[20]. For further insight on the alternative methods in calculating the orbital precession and planetary orbits in Schwarzschild, Reissner-Nordstrom, Kerr and Kerr-Newman space-times see [21]–[29].

In this paper, instead of using Einstein’s direct method and writing the integrand in term of the small parameter $\frac{GM}{r}$ or using Newtonian and Post-Newtonain approximation, we consider a particulary simple geodesic, namely the circular orbit. We then add the first, second and higher-order deviations to such an orbit, thereby obtaining the exact trajectory of the neutral test particles. The advantage of this approach is in finding explicit form of the trajectory for arbitrary values of parameter $\frac{GM}{r}$ i.e., in the vicinity of very massive and compact bodies, and also to take into account the full relativistic effects without perturbing the curvature of the space-time. However, this method is valid for small eccentricities. In addition to the method presented here, generalization of geodesic deviation equations to second and third orders are presented in [30]-[32]. Also in [33], a derivation based on covariant bitensor formalism [34]-[36] for higher orders geodesic deviation has been presented.

The paper is arranged as follows. In section 2 we give a brief review of the method used here. The circular motion in Kerr-Newman metric is explained in section 3. In section 4, we study the first geodesic deviation around this circular orbit. To obtain the trajectory of test particles in Kerr-Newman field, we derive the second and third-order geodesic deviations in section 5. Finally, conclusions are presented in section 6.

2 Geodesic deviation method

One of the most important equations in General Relativity is the geodesic deviation equation which, for a neutral particle is given by [7, 37]

[TABLE]

As is shown in the Figure 2, $u^{\alpha}(\tau,\lambda)=\frac{\partial x^{\alpha}}{\partial\tau}$ is a time-like tangent vector along the geodesic and $n^{\alpha}(\tau,\lambda)=\frac{\partial x^{\alpha}}{\partial\lambda}$ is the deviation four-vector. Also, $\frac{D}{D\tau}$ is the covariant derivative along the proper time $\tau$ and $R^{\mu}_{\,\,\,\,\lambda\nu\kappa}$ is the Riemann tensor of the space-time. It is often convenient to rewrite the above equation in the non-trivial covariant form

[TABLE]

Now, we use the geodesic deviation method to construct one geodesic from an initial one by expansion in a series of powers of small parameter [8]

[TABLE]

Our aim is to obtain an expression in terms of deviation vectors. In the above equation, the second term, $\frac{\partial x^{\mu}}{\partial\lambda}$, is the definition of deviation vector and shows deviation from first-order geodesic. However, in third term, $\frac{\partial^{2}x^{\mu}}{\partial\lambda^{2}}$, is not a vector anymore. One is therefore required to change this term to an expression showing the second-order geodesic deviation. To this end, we define the vector $b^{\mu}$

[TABLE]

Now, by substituting this in equation (5), one obtains an expression in terms of the order of vector deviations

[TABLE]

where $\epsilon=(\lambda-\lambda_{0})$. For simplicity we introduce $\delta^{n}x^{\mu}(\tau)$ as the $n$-th order of deviation, then equation (7) takes the form

[TABLE]

where $\delta x^{\mu}(\tau)=\epsilon n_{0}^{\mu}(s)$ and $\delta^{2}x^{\mu}(\tau)=\epsilon^{2}(b^{\mu}-\Gamma^{\mu}_{\nu\lambda}n^{\nu}n^{\lambda})$ are first and second-order geodesic deviations, respectively. With straightforward calculations, one finds that vector $b^{\mu}$ satisfies

[TABLE]

Rewriting the above equation in a non-trivial covariant form, we have

[TABLE]

which is called the second-order geodesic deviation equation. Generalization of this equation for the case of the charged particles has recently been obtained in [38]. It is worth noting that the left-hand sides of equations (4) and (10) are the same, but in the latter there are some extra terms. It can easily be shown that by extending the method to third and higher orders, this general pattern persists. In other words, the geodesic deviation equations of higher order have the same left-hand side and differ from each other in the extra terms appearing on the right-hand side which depend on the solutions of lower orders.

In the next sections, by considering a circular orbit of radius $R$ as the original geodesic $x^{\mu}_{0}(\tau,\lambda)$ for a neutral particle in Kerr-Newman space-time and using equations (4) and (10), we obtain the components of deviation vectors $n^{\mu}$ and $b^{\mu}$ in such space-times.

3 Circular orbits in Kerr-Newman metric

Let us now study the circular orbits outside a massive, charged, rotating source [40]. This space-time which is described by the Kerr-Newman metric is an exact, axis-symmetric and asymptotically flat solution of the Einstein-Maxwell equations. In the Boyer$-$Lindquist coordinates (with $c=G=1$) the line-element is given by

[TABLE]

with

[TABLE]

where parameters $M$, $a$ and $Q$ stand for mass, angular momentum per unit mass and charge of the black-hole. As is clear, in the limiting cases $Q=0$ and $a=0$ the Kerr-Newman metric reduces to the Kerr and Reissner-Nordstrom metrics. Also it reduces to the Schwarzschild metric for $a=0$, $Q=0$. For further insight, see [37].

Since the Kerr-Newman metric has no explicit dependence on $t$ and $\varphi$ coordinates, two conserved quantities, namely energy ($\varepsilon$) and angular momentum ($l$) can be obtained from Euler-Lagrange equations

[TABLE]

[TABLE]

where a dot represents derivative with respect to the proper time $\tau$. For the radial coordinate $r$ we find

[TABLE]

from which we obtain a relation for the radial acceleration

[TABLE]

Now, let us consider a neutral test particle with mass $m\ll M$ along a circular orbit of radius $R$ in the equatorial plane. Then using constraints $\theta=\frac{\pi}{2}$ and $\dot{r}=\ddot{r}=0$ in the above equations, components of the four-velocity vector are given by

[TABLE]

It is clear that with increasing charge and rotation parameters the angular frequency of the circular orbit (unperturbed initial orbit), $\omega_{0}$, decreases.

4 The first-order geodesic deviation around circular orbits

Now, we are ready to solve the first-order geodesic deviations in a Kerr-Newman space-time. To this end, in equation (4) the radial coordinate $r$ is set to constant $R$, the four-velocities are given by equation (17), and our calculations are restricted to the equatorial plane, $\theta=\frac{\pi}{2}$. Thus the first-order geodesic deviation equation takes the form

[TABLE]

where

[TABLE]

with

[TABLE]

From equation (18), we can see that $n^{\theta}$ component of geodesic deviation equation is independent of variables $n^{t}$,$n^{r}$ and $n^{\varphi}$. The harmonic oscillator equation for $n^{\theta}$ has an angular frequency given by

[TABLE]

with the following solution

[TABLE]

In the Schwarzschild and the Kerr limit, the above result reduces to that of references [8] and [11], respectively. It is worth noting that in the Schwarzschild limit we have $\omega_{\theta}=\omega_{0}$, namely the period is equal to the period of the planetary motion which can be interpreted as a pure coordinate effect so in this case we can take $n^{\theta}=0$.

Now, by eliminating derivatives of $n^{\varphi}$ and $n^{t}$ in differential equation for $n^{r}$, we find the following harmonic oscillator equation

[TABLE]

with the characteristic frequency

[TABLE]

where

[TABLE]

[TABLE]

resulting in perihelion and aphelion corresponding to $\tau=0$ and $\tau=\frac{\pi}{\omega}$, respectively. We note that the discrepancy between the two angular frequencies $\omega_{0}$ and $\omega$ is what causes the perihelion precession. The solutions $n^{\varphi}$ and $n^{t}$ are given by

[TABLE]

where the amplitudes depend on $n_{0}^{r}$

[TABLE]

Now, these first-order deviations (linear approximation) should be added to the unperturbed circular orbit (zeroth-order approximation) to obtain the new trajectory as follows

[TABLE]

[TABLE]

[TABLE]

[TABLE]

If we compare equation (29) with the Kepler result [41]

[TABLE]

we may identify the eccentricity $e$ with $\frac{n_{0}^{r}}{R}$ and the greater half-axis $a$ with $R$. Also, we can obtain an approximation for the orbital precession of neutral test particles in the equatorial plane of the Kerr-Newman space-time, using first-order deviations

[TABLE]

As one might expect, for the most general black-hole, i.e. Kerr-Newman, in addition to the mass, the precession of perihelion is also modified by the angular momentum and electric charge of the source. In the method presented here, the orbital precession contains higher-order terms, $\frac{M}{R}$, which demonstrate higher accuracy as compared to that of Einstein’s direct method. In the lowest-order of approximation the orbital precession can thus be written as

[TABLE]

It is now interesting to investigate how different parameters of a central body influence the perihelion advance. The first term is the standard value for the perihelion precession of a test particle in the Schwarzschild field. The second and third terms are the effect of the rotation of central body. If the rotation of central body and the orbital motion of test particle are in the same or opposite directions, the linear term in the rotation parameter $a$ can make the perihelion advance to increase or decrease, respectively. Also, the quadratic term in $a$ always leads to higher precession since such extra energy would add extra mass to the central body according to the equivalence principle. Moreover, an extra term for the Reissner-Nordstorm space-time is observed and since the charge parameter appears in the quadratic terms, the sign of $Q$ is not important and the charge of the black-hole causes the precession to diminish (the charge parameter of the black hole seems to induce repulsive effects [42]). As can be seen from the above equation, orbital precession contains a new term $\frac{216\pi a\sqrt{M}Q^{2}}{R^{7/2}}$, namely the Kerr-Newman effect, which couples the gravitational mass of the source, angular momentum and charge and tends to enhance the perihelion advance.

In the absence of charge, equation (34) reduces to the usual orbital precession in Kerr space-time [11]. Also if the rotation of the gravitational source vanishes and the source becomes neutral, the result reduces to that of the perihelion precession in a Schwarzschild field, that is

[TABLE]

which corresponds to equation (48) of reference [8]. By comparing the above equation to equation (2), we find extra quadratic and cubic terms of the form $\frac{M}{R}$ which would improve the accuracy. It also shows that the geodesic deviation method can be used in the vicinity of massive and compact objects.

It is also clear that to first-order deviations we only keep terms up to $e^{2}$. In what follows we derive and solve higher-order geodesic deviations and obtain non-linear terms involving $e^{2}=(\frac{n_{0}^{r}}{R})^{2}$, together with other higher-order effects.

5 Second and third-order geodesic deviation

For calculating the trajectory of neutral test particles in a Kerr-Newman space-time, we need to obtain a set of deviation vectors $(n^{\mu}(\tau),b^{\mu}(\tau),\cdots)$ on the reference geodesic. To this end, we insert the complete solution of the first-order deviation vector, equations (22), (26), and (27) into equation (10), leading to a set of linear second-order differential equations for deviation vector $b^{\mu}(\tau)$

[TABLE]

where the coefficient $\epsilon$ shows dependence of $\delta^{2}x^{\mu}$ on the amplitude $n_{0}^{r}$ ($n_{0}^{\theta}$), and constants $C^{r}$, $C^{\theta}$, $C^{\varphi}$ and $C^{t}$ contain quantities depending on $M$, $R$, $a$, $Q$ and $\omega$ , $\omega_{0}$ and periodic functions of the form $\sin(2\omega\tau)$, $\cos(2\omega\tau)$ $\sin(2\omega_{\theta}\tau)$ and $\cos(2\omega_{\theta}\tau)$

[TABLE]

[TABLE]

[TABLE]

[TABLE]

The solution for $b^{\mu}(\tau)$ has two terms, the general oscillatory solution with angular frequency $\omega$ and $\omega_{\theta}$ that accounted for $n^{\mu}$, and the second particular solution which includes the oscillatory terms with angular frequency $2\omega$, $2\omega_{\theta}$, $(\omega+\omega_{\theta})$ and $(\omega-\omega_{\theta})$ as follows

[TABLE]

[TABLE]

[TABLE]

[TABLE]

To obtain the geodesic curve $x^{\mu}$, we calculate $\frac{1}{2}\delta^{2}x^{\mu}$ with the result

[TABLE]

[TABLE]

[TABLE]

[TABLE]

Coefficients in the above equations can be obtained after long calculations and are listed in the Appendix. As we expect in the Schwarzschild limit $a=0$, $Q=0$, and Kerr limit $Q=0$, the solution reduces to that of references [8] and [11], respectively. By adding the second-order deviation vector $b^{\mu}(\tau)$ (with angular frequency $2\omega$) to the circular orbit we can get a better approximation of the elliptical orbit. The resulting trajectory is not an ellipse since $R\neq a$ and $e\neq\frac{n_{0}^{r}}{R}$, but for $x^{\mu}$ the perihelion and aphelion can be derived from $\omega\tau=2k\pi$ and $w\tau=(1+2k\pi)$ where $k\in Z$, respectively. Therefore the semi-major axis $a$ and eccentricity $e$ of an ellipse are given by

[TABLE]

[TABLE]

Finally, $\tau$ can be substituted in $r(\tau)$, resulting in $r(\varphi)$ up to the second order

[TABLE]

In the Schwarzschild limit, we have an elliptical orbit

[TABLE]

which is different from analogous result [8] in a constant term, because we have used another initial conditions. Then for the above result we obtain

[TABLE]

and

[TABLE]

Since the variable $u(\varphi)=\frac{1}{r}$ is more convenient than $r$, we compare the numerical solution of the relativistic version of Binet’s equation and its approximate solution for $u(\varphi)$ with the $u(\varphi)$ obtained via the higher-order geodesic deviation method, equation (53), in figure 3. As can be seen, our result has more coincidence with numerical graph.

Also for the Schwarzschild case the shape of the orbit in the limit $\frac{M}{R}\rightarrow 0$ is described up to second-order in $e=\frac{n_{0}^{r}}{R}$ by

[TABLE]

which is an exact equation of an ellipse with $a\simeq R(1-\frac{e^{2}}{2})$, in agreement with equations (50) and (51). In the third-order approximation, one may find a covariant deviation equation for $\delta^{3}x^{\mu}$ (see equation (72) of [8]), which can be written in matrix form. As previously mentioned, this matrix has the same left-hand side as those of equations (18) and (37), but a different vector in the right-hand side. Components of this vector contain terms including sine and cosine functions of the form $\omega\tau$ ($\omega_{\theta}\tau$) and $3\omega\tau$ ($3\omega_{\theta}\tau$). The presence of the functions $\cos(\omega\tau)$ and $\cos(\omega_{\theta}\tau)$, are the same as resonance terms, leading to solutions proportional to $\tau$. Since such terms in the solutions correspond to unbounded deviations, we can use Poincare’s method to avoid this problem [8].

6 Conclusions

In this paper we have generalized the geodesic deviation method [8] to the Kerr-Newman space-time which represents the most general black-hole. Our method is fully relativistic and constitutes an alternative approach for calculating the perihelion precession and trajectory of test particles in gravitational fields. We started with a circular orbit with constant velocities. In this case the first-order geodesic deviations, equations (4), was reduced to a system of linear differential equations with constant coefficients which is a collection of harmonic oscillators, equation (18). The discrepancy between two circular frequencies, the unperturbed circular orbit $\omega_{0}$ and the perturbed orbit $\omega$ to first order, equations (17) and (24), is exactly what produces the precession of planet’s perihelion. The effect of mass of a black-hole, the angular momentum and electric charge of the source on the orbital precession was also discussed. We then proceeded to derive the second-order geodesic deviations, equation (10), leading to a collection of harmonic oscillators excited by the external periodic forces represented by the right hand side of equation (37), which may be solved easily. Finally, by adding the first corrections (linear approximation), second and higher-order corrections (non-linear approximation) to the circular orbit (zeroth-order approximation), we constructed the trajectory of neutral test particles, equation (52), in the Kerr-Newman space-time. Also, the study of non-equatorial orbits and the precession of the orbital plane in such a space time, using geodesic deviation method, would be an extension of our results for a future work.

Finally, it is worth stressing that in contrast to the Post-Newtonian expansion [44, 45, 46] where the orbits are calculated in an almost Minkowski space-time and the corrections due to the non-zero curvature are added perturbatively, there is no such constraints on the curvature of space-time in the geodesic deviation method and the fully relativistic effects are taken into account. Also, the additional terms of the form $M/R$ obtained in the geodesic deviation method when compared to similar solutions in the standard approximation result in better accuracies. The other advantage of the geodesic deviation method is to provide a way for calculating small-eccentricity orbits in the vicinity of massive and compact objects. This is so because in the Post-Newtonian approximation, the assumption of weak gravitational fields and small velocities would not be satisfied in scenarios where the ratio $M/R$ is large.

Acknowledgements

We would like to thank Richard Kerner for helpful discussions and Mohsen Khodadi for a careful reading of the manuscript and useful suggestions. We also thank the anonymous referee for valuable comments.

Appendix

The coefficients of $\delta^{2}x^{\mu}$ in the second-order geodesic deviation equations are

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]