Direct and Exact Description of Null Geodesics in Schwarzschild Spacetime

Leo Morgovsky

Johann Wolfgang Goethe-Universität Frankfurt am Main, Frankfurt, Germany.

DOI: 10.4236/jhepgc.2023.91020                

Abstract

Null geodesics for massless particles in Schwarzschild spacetime are obtained by direct integration of the trajectory equation in spatial coordinates without transformation to the inverse space. The results are classified following Chandrasekhar depending on the ratio of the impact parameter of the trajectory to its critical value. In the subcritical and supercritical cases the geodesics are expressed in terms of elliptic integrals of the first kind. Some results are formally different from the classical ones, but in fact equivalent to them, being at the same time more compact and descriptive.

Keywords

Black Holes, Schwarzschild Metric, Photon Trajectory, Hawking Radiation

1. Introduction

Despite the fact that the study of recently discovered gravitational waves offers a new method of obtaining knowledge about the black holes (BH), the analysis of the trajectories of massive and massless particles is still the most important source of information about these largely mysterious physical objects.

Recently obtained by the collaboration of the event horizon telescope (EHT) data on light bands in the vicinity of the supermassive object in the heart of the galaxy M87 [1] once again demonstrated the importance of solutions for null geodesics in the Schwarzschild and Kerr spacetimes. Despite the fact that these solutions have been obtained for a long time, efforts continue to find alternative approaches that would make it easier to obtain or clearer the final expressions. The most significant pioneering works in the field have been carried out by C. Darwin [2] [3] (in Schwarzschild metric) and S. Chandrasekhar [4] (in Kerr metric). The results of their works are described in detail in Chandrasekhar’s own book as well as in books of other authors [5] [6]. A review over the last decade can be found in Ref. [7] ).

All calculations accomplished so far have been performed through transformation of radial coordinate r into an auxiliary intermediate function u = 1/r. This technique originates in classical celestial mechanics long before the relativistic era and is explained by quite clear reason. It is a desire to remove the singularity at the origin of coordinates to infinity and, at the same time, to include the physical infinity, which is a source and (often) an outlet for test bodies, and very often the position of the observer, in the field of events.

Most often, researchers chose the use of Jacobi elliptic functions [8] [9] [10] as a solution method. More recently, it has been shown that Weierstrass elliptic function can also be successfully used in these calculations [11] [12]. The advantage of both approaches is the opportunity to get results in a generally accepted form r = r(φ).

The solution in the form of elliptic integrals is obtained as an inverse form, which was often considered a disadvantage. However, when calculating the trajectory on a computer, it does not matter which of the forms is the preferred one. Both are completely equivalent and can easily be converted into one another graphically and/or numerically.

The aim of this work was to solve the equation of null geodesics using the natural coordinate system radius-azimuth angle, without additional coordinate transformations. In this case, the solution seems to be more natural and direct. This made it possible to obtain solutions that sometimes formally differ from the well-known ones, but in fact completely equivalent to them and, at the same time, are often simpler and more illustrative.

This approach has an advantage due to its simplicity and straightforwardness in presentation.

The consideration of the problem of searching for null geodesics in this paper is focused exclusively on the most promising geometry for further applications, when the source of test particles is localized at a sufficiently large distance from the event horizon, at least much greater than the Schwarzschild radius.

2. Exact Solutions in Terms of Elliptic Integrals

The Schwarzschild metric element has the form (c = 1),

$\text{d}{s}^2=-f\left(\rho \right)\text{d}{t}^2+f{\left(\rho \right)}^{-1}\text{d}{\rho }^2+{\rho }^2\left(\text{d}{\theta }^2+{\sin }^2\theta \text{d}{\phi }^2\right)$, (1)

where $f\left(\rho \right)=1-R_s/\rho$, $R_s=2GM$ the Schwarzschild radius, G is the gravitational constant and M is the mass of the body which creates gravitational field.

First of all, the solution of the problem is simplified by the fact that for the massless particles moving at the speed of light, the square of the interval $\text{d}{s}^2=0$. In addition, spherical symmetry of the Schwarzschild field makes it possible to consider the plane of particle motion as equatorial one ( $\theta =\pi /2$ ). Finally, it should be taken into account the conservation laws of angular momentum L and energy E. As a result of all these circumstances from the metric (1) follows the equation [9] [10]

${\left(\frac{\text{d}r}{\text{d}\phi }\right)}^2=\frac{r^4}{b^2}-r^2+r$, (2)

where we have introduced a dimensionless coordinates $r={\rho /R}_s$ and an impact parameter b = L/ER_s.

Integral of this equation will allow us to determine the required trajectories

$\phi =\pm b{\displaystyle \int \frac{\text{d}r}{\sqrt{r{P}_3\left(r\right)}}}$, (3)

where P₃(r) is a cubic polynomial.

P₃ has 3 roots with the property that their sum ${\displaystyle \underset{i}{\sum }{r}_i}=0$. This property is a

direct consequence of the absence of a quadratic term in the polynomial. Depending on the sign of the discriminant of the cubic polynomial

$D=\frac{b^4}{4}\left(1-\frac{4b^2}{27}\right)$ (4)

there are three possible types of roots which determine the physical properties of null geodesics. All of them differ depending on the value of the impact parameter regarding to its critical value at which the discriminant (4) turns to 0. It

follows from Equation (4) that the critical impact parameter $b_c=\frac{3\sqrt{3}}{2}$.

For b = b_c, the polynomial has 3 real roots, two of which are equal, for b > b_c all roots are real and different, and finally, for b < b_c, only one root is real, and two others are complex conjugate.

We will consider all these cases separately.

1) b = b_c.

In this case two roots are equal and the equation becomes

$\phi =\pm b_c{\displaystyle \int \frac{\text{d}r}{\left(r-r_1\right)\sqrt{r\left(r-r_3\right)}}}$. (5)

The calculation of this integral is given in the Appendix A.

As a result

$r\left(\phi \right)=\frac{3}{2}+\frac{9}{4\left[\cosh \left(\phi \right)-1-\frac{\sqrt{3}}{2}\sinh \left|\phi \right|\right]}$. (6)

It can be seen from the Equation (6) that the (unstable) critical photon trajectory is circular and its radius asymptotically approaches one and a half Schwarzschild radius (in a system of units, where G = 1, the result for the asymptotic value has the form r = 3M). Approximation to the asymptote occurs exponentially, approaching with every full turn the radius 3/2 as exp(−2π).

Qualitatively, these results are well known and discussed many times. However the solution found by C. Darwin [2] [3] in terms of reversal coordinate u = 1/r (discussed in details by S. Chandrasekhar in his classical book [4] ) has the form

$u=-\frac{1}{3M}+\frac{1}{2M}{\tanh }^2\left[\left(\phi -{\phi }_0\right)/2\right]$, (7)

where integration constant ${\tanh }^2\left({\phi }_0/2\right)=1/3$.

Both solutions have the same properties and are obviously different. It can be shown nevertheless by algebraic transformations that these solutions are identical. A more detailed analysis is presented in the work of the author [13].

Some particular interest is the case when a photon with critical impact parameter does not appear from infinity, but is born near the event horizon, for example, as a particle of Hawking radiation [6]. In this case the Equation (5) remains valid, but boundary condition now has the form r(0) = 1.

A solution completely analogous to that described in Appendix A leads to the following result

$r=\frac{3}{2}-\frac{9}{14\cosh \left(\phi \right)+8\sqrt{3}\sinh \left|\phi \right|+4}$, (8)

Thus, the photon tends to the same critical orbit 3/2, but from the inner side of the spherical surface and also approaches it, spinning in a spiral.

2) b< b_c.

In this case one root is real while two others are complex conjugate.

The geodesic trajectory now defined by the integral

$\phi =\pm b{\displaystyle \underset{r}{\overset{\infty }{\int }}\frac{\text{d}r}{\sqrt{r\left(r-r_1\right)\left(r-r_2\right)\left(r-r_3\right)}}}$, (9)

where the integration limits are chosen in such a way as to automatically satisfy the boundary condition φ = 0 as r → ∞.

After transformation $\left(r-r_2\right)\left(r-r_3\right)=\left[{\left(r-\mathrm{Re}{r}_2\right)}^2+\mathrm{Im}{r}_2^2\right]$ the result of integration can be formally written according to Ref. [14] (Eq. 260.00) as expressed in terms of

$\phi =\pm \frac{b}{\sqrt{AB}}\left[F\left({\psi }_{\infty },k\right)-F\left(\psi ,k\right)\right]$. (10)

where F(ψ, k) is an incomplete elliptic integrals of the 1^st kind.

In addition, we have introduced the following notation

$A=\left|r_2\right|$, $B=\sqrt{{\left(r_1-\mathrm{Re}\left(r_2\right)\right)}^2+{\mathrm{Im}}^2\left(r_2\right)}$, $k^2=\frac{{\left(A+B\right)}^2-r_1^2}{4AB}$ (11)

$\psi =\arccos \left[\frac{\left(A-B\right)r-r_{1}A}{\left(A+B\right)r-r_{1}A}\right]$, ${\psi }_{\infty }=\arccos \left[\frac{A-B}{A+B}\right]$ (12)

where the roots are defined with expressions

$r_1=-\frac{2b}{\sqrt{3}\sin \left(2\alpha \right)}$, $\mathrm{Re}\left(r_2\right)=\frac{b}{\sqrt{3}\sin \left(2\alpha \right)}$, $\mathrm{Im}\left(r_2\right)=b\cot \left(2\alpha \right)$, (13)

$\alpha =\arctan \sqrt[3]{\tan \left(\beta /2\right)}$, $\sin \beta =b/b_c$. (14)

Trajectories of the particles whose impact parameter is less than the critical one inevitably terminate inside the event horizon. The result (10) obtained here differs from those found by other authors, since no one required fullfilment of natural condition φ = 0 at infinity.

The corresponding trajectories are depicted in (Figure 1) for various b < b_c (1 and1a b = 2.55; 2 and 2a b = 2; 3 and 3a b = 1; 4 and 4a b = 0.5). Trajectories 1-4 belong to positive φ while symmetric to them 1a-4a belong to φ < 0.

3) b> b_c.

All three roots are real and different in this particular case. The integral that determines the geodesic trajectory has the same form as Equation (9),

$\phi =\pm b{\displaystyle \underset{r_1}{\overset{r}{\int }}\frac{\text{d}r}{\sqrt{r\left(r-r_1\right)\left(r-r_2\right)\left(r-r_3\right)}}}$, (15)

in which, however, the lower limit of integration is the coordinate of a periapsis point and upper limit is a particle coordinate. The equation which defines the first one is the same trajectory Equation (2) where is the maximum point of the trajectory, i.e. the point of the nearest approach of the particle to the gravitating body. It is quite clear that at this point $\text{d}r/\text{d}\phi =0$. The roots in this case are

$r_1=\frac{2b}{\sqrt{3}}\cos \frac{\alpha }{3}$ ; $r_{2,3}=\frac{2b}{\sqrt{3}}\cos \left(\frac{\alpha \pm 2\pi }{3}\right)$ ; $\alpha =\pi -\arccos \left(\frac{b_c}{b}\right)$. (16)

Accordingly, the periapsis coordinate is the largest root $r_1$.

The result of integration is [14] (Eq.258.00)

$\phi \left(r\right)=\pm \frac{\sqrt{2}\sqrt[4]{3}}{\sqrt{\sin \frac{2\alpha }{3}}}F\left(\psi ,k\right)$, (17)

where $k^2=\frac{2\sin \left(\frac{\pi +\alpha }{3}\right)\cos \left(\frac{2\pi -\alpha }{3}\right)}{\sin \left(2\alpha /3\right)}$, and (18)

Figure 1. Particle trajectories with b < b_cdestined to fall inside the event horizon. 1, 1a − b = 2.55; 2, 2a − b = 2; 3, 3a − b = 1; 4, 4a − b = 0.5; 1 − 4 φ > 0; 1a − 4a φ < 0.

${\sin }^2\psi \left(r\right)=\frac{\left[r-\frac{2b}{\sqrt{3}}\cos \left(\alpha /3\right)\right]\sin \left(\alpha /3\right)}{\left[r-\frac{2b}{\sqrt{3}}\cos \left(\frac{2\pi -\alpha }{3}\right)\right]\sin \left(\frac{\pi +\alpha }{3}\right)}$. (19)

The corresponding trajectories are shown in Figure 2. The particles fly around the black hole and then return back to infinity. Due to the symmetry with respect to the horizontal plane, trajectories +φ and −φ are actually the same ones. The dashed line shows the position of the critical photon orbit.

As in previous cases, the resulting Equation (17) formally differs from the well-known classical solution [4]. Despite this, the calculations show complete and exact equivalence of both solutions. The proof of the identity of both solutions is given in Appendix B. The advantage of the above obtained solution is its complete closure. At the same time, the classical solution [4] uses the periapsis value, which must be introduced from outside the equation, for example, from astronomical data. No additional data besides the impact parameters is necessary in our approach.

These results make it possible to write out a closed equation for determining the photons deviation angle after turning around the BH.

It is quite obvious that the angle of total deviation $\Omega$ is expressed in terms of the deflection angles as

$\Omega =2\left|{\phi }_{\infty }\right|-\pi$, (20)

where ${\phi }_{\infty }$ is a deflection angle at infinity.

Using the expression (17), the desired result takes the form

$\Omega =\frac{2\sqrt{2}\sqrt[4]{3}}{\sqrt{\sin \frac{2\alpha }{3}}}F\left({\psi }_{\infty },k\right)-\pi$, (21)

where the angular variable in the integral is found from the formula (19) and has

Figure 2. Particle trajectories withb > b_cwhich are deflected by the BH. 1—b = 1.1b_c; 2—b = 1.25b_c; 3—b = 1.5b_c; 4—b = 2.0b_c; 5—critical photon orbit.

a simple form

${\sin }^2{\psi }_{\infty }=\frac{\sin \alpha /3}{\sin \frac{\pi +\alpha }{3}}$. (22)

The Equation (21) makes time-consuming series expansion $\Omega \left(b\right)$ [10] unnecessary. The hyperbolic curve of this dependence is shown in Figure 3. The equation which expresses this dependence can be useful in the analysis of the effects of gravitational lensing, one of the most frequently discussed topics.

The listed options exhaust all possible cases of trajectories of massless particles starting at infinity. It is useful and instructive, however, to consider several special cases on the basis of the obtained solutions.

Consider the case of a very weak gravitational field, in which $b\gg b_c$. This situation occurs in particular when the photon beams graze the edge of the solar disk. In this case we conclude from Equation (16) that $\alpha \approx \pi /2$. Expanding Equations (18) and (22) to the terms of the first non-vanishing order, we obtain

$k^2=2/b$ and ${\psi }_{\infty }=\frac{\pi }{4}+\frac{3}{4}{b}^{-1}$.

Expanding then Equation (17), including the integral and taking into account that $k^2\ll 1$, we get ${\phi }_{\infty }=\frac{\pi }{2}+b^{-1}$. From Equation (21) we finally obtain $\Omega =\frac{2}{b}$ (or in dimensional units it is twice the ratio of gravitational and real physical radii of the Sun 2Rs/Rsol).

This result was first obtained personally by A. Einstein and subsequently confirmed by the expedition of the Royal Astronomical Society during the solar eclipse. As it was shown in Ref. [10], this is only the first term in the expansion of the deviation angle in a power series of the impact parameter. However, for the relatively weak gravitational field of the Sun, this turned out to be enough.

Figure 3. Plot of the deviation angle depending on the impact parameter.

It is clearly seen from Figure 3, as the impact parameter decreases to the critical value, the deviation angle increases and for an exponentially small difference between them, a photon is able to make several turns around the BH. The closer the photon’s orbit approaches the critical circular orbit, the more rotations the photon can make. Calculation results according to Equations ((21), (22)) are shown in Figure 4. The number of rotations N = Ω/2π is expressed as a function of the exponent in the difference between the impact parameter and its critical value log(b/b_c−1).

Figure 4. Number of rotations depending on the difference between the impact parameter and its critical value.

Figure 5. Photon trajectory computed for b/b_c − 1 = 10⁻³ (The distance scale is in the units of Schwarzschild radius).

The trajectory calculated by approximate formula according to Eqs. (17)-(19) for b/b_c = 1 + 10⁻³ is shown in Figure 5.

Appendix A

In this case, the roots of the cubic polynomial are easily found by elementary formulas and they are equal $r_1=r_2=\frac{3}{2}$ while $r_3=-3$.

The integral (5) is calculated using the substitution $t={\left(r-\frac{3}{2}\right)}^{-1}$. As a result it takes the form

$-{\displaystyle \int \frac{\text{d}t}{\sqrt{1+6t+\frac{27}{4}{t}^2}}}$. (A.1)

This integral is calculated in a standard way and is equal to

$-\frac{2}{\sqrt{27}}\ln \left(\frac{\sqrt{27\left(1+6t+\frac{27}{4}{t}^2\right)}+\frac{27}{2}t+6}{C}\right)$, (A.2)

where C is a constant of integration.

Substituting this result into the formula (5) we get the equality

$\phi =-\ln \left(\frac{\sqrt{\frac{27}{4}\left(1+6t+\frac{27}{4}{t}^2\right)}+\frac{27}{2}t+6}{C}\right)$. (A.3)

Natural boundary condition is that, when r → ∞ that means t = 0 φ = 0.

Hence it follows that $C=6+3\sqrt{3}$.Substituting C and expressing now t as a function of φ, we obtain

$t=\frac{\left(2+\sqrt{3}\right)\exp \left(-\phi \right)+{\left(2+\sqrt{3}\right)}^{-1}\exp \left(\phi \right)-4}{9}$. (A.4)

Returning back to the variable r, we finally find

$r=\frac{3}{2}+\frac{9}{4\left(\cosh \left(\phi \right)-1-\frac{\sqrt{3}}{2}\sinh \left|\phi \right|\right)}$. (A.5)

Appendix B

We will outline the course of the proof for the case b > b_c.

The well-known formula [4] for supercritical geodesics has the form

$\phi =2\sqrt{\frac{P}{Q}}\left[K\left(k\right)-F\left(\frac{\chi }{2},k\right)\right]$, (B.1)

where P is a periapsis and $Q=\sqrt{\left(P-2M\right)\left(P+6M\right)}$.

It is quite easy to prove that in our notation they are equal to $P=r_1$ and $Q=\sqrt{\left(r_1-1\right)\left(r_1+3\right)}$ because $\text{2}M=R_s=\text{1}$ in our units system.

After that it is easy to verify that the factor in front of the integrals in Equation (A2.1) is exactly the same as the factor in front of the integral in the Equation (17). The modules of the integrals k also coincide. Finally, the equality of the integrals

$K\left(k\right)-F\left(\frac{\chi }{2},k\right)=F\left(\psi ,k\right)$, (B.2)

turns out to be an identity valid for any values of variables r and b.

Conflicts of Interest

The author declares no conflicts of interest regarding the publication of this paper.

References

[1]	Johnson, M.D., Lupsasca, A., Strominger, A., et al. (2020) Universal Interferometric Signatures of a Black Hole’s Photon Ring. Science Advances, 6, Article No. eaaz1310. https://doi.org/10.1126/sciadv.aaz1310
[2]	Darwin, C. (1959) The Gravity Field of a Particle I. Proceedings of the Royal Society A, 249, 180-194. https://doi.org/10.1098/rspa.1959.0015
[3]	Darwin, C. (1961) The Gravity Field of a Particle II. Proceedings of the Royal Society A, 263, 39-50. https://doi.org/10.1098/rspa.1961.0142
[4]	Chandrasekhar, S. (1992) The Mathematical Theory of Black Holes. Clarendon Press, Oxford.
[5]	Misner, C.W., Thorn, K.S. and Weeler, J.A. (1973) Gravitation. Freeman, San Francisco.
[6]	Shapiro, S.L. and Teukolsky, S.A. (1983) Black Holes, White Dwarfs and Neutron Stars. The Physics of Compact Objects. Whiley-Interscience Publication, New York. https://doi.org/10.1002/9783527617661
[7]	Hackmann, E. and Lämmerzahl, C. (2014) Analytical Solution Methods for Geodesic Motion. AIP Conference Proceedings, 1577, Article No. 78. https://doi.org/10.1063/1.4861945
[8]	Luminet, J.-P. (1979) Image of a Spherical Blake Hole, Astron. Astrophys, 75, 228-235.
[9]	Cadez, A. and Kostic, U. (2005) Optics in the Schwarzschild Space-Time. Physical Review D, 72, Article ID: 104024. https://doi.org/10.1103/PhysRevD.72.104024
[10]	Munoz, G. (2014) Orbits of Massless Particles in the Schwarzschild Metric: Exact Solutions. American Journal of Physics, 82, 564-573. https://doi.org/10.1119/1.4866274
[11]	Scharf, G. (2011) Schwarzschild Geodesics in Terms of Elliptic Functions and Related Red Shift. Journal of Modern Physics, 2, 1-21.
[12]	Gibbons, C.W. and Vyska, M. (2012) The Application of Weierstrass Elliptic Functions to Schwarzschild Null Geodesics. Classical and Quantum Gravity, 29, Article ID: 065016. https://doi.org/10.1088/0264-9381/29/6/065016
[13]	Morgovsky, L. (2022) Critical Equatorial Trajectories of Massless Particles in a Kerr spacetime. Journal of High Energy Physics, Gravitation and Cosmology, 8, 647-656. https://doi.org/10.4236/jhepgc.2022.83047
[14]	Bird, P.F. and Friedman, M.D. (1971) Handbook of Elliptic Integrals for Engineers and Scientists. 2nd Edition, Springer-Verlag, New York. https://doi.org/10.1007/978-3-642-65138-0