# Light ray paths near schwarzschild blackhole

1. Oct 22, 2009

### nocks

hey there, i'm interested in (eventually) simulating light ray paths near black holes, starting with schwarzschild blackholes and working my way to kerr-newman blackholes.
I have a good understanding of the nature of blackholes but have trouble when it comes to the equations.
My background is in computer science and I was wondering if anyone here could put a non-physicist on the right track to plotting the trajectory of a light ray near a black hole. ( in 2d for now ).
Thanks

Last edited: Oct 23, 2009
2. Oct 24, 2009

### nocks

I'm currently going over lagrangian mechanics

3. Oct 24, 2009

### ZikZak

Last edited by a moderator: Apr 24, 2017
4. Oct 24, 2009

### Naty1

5. Oct 24, 2009

### nocks

Thanks for the paper, all the others I have read from arxiv have not been very non-physicist friendly.

Last edited by a moderator: Apr 24, 2017
6. Oct 25, 2009

### George Jones

Staff Emeritus
Are you interested in what an observer would see with his eyes, or in plotting the trajectory of a photon on a coordinate map, or maybe both? The first thing is somewhat hard, and the second thing is easier.

Let me use an example to illustrate what I mean. Suppose a ship makes a long ocean voyage. The voyage could be watched from a telescope on a satellite in geosynchronous orbit, or the voyage could be plotted as a moving, glowing dot on a map.

7. Oct 25, 2009

### nocks

My end goal will be to render what an observer would see , but for now I would like to simply plot the trajectory of several photons on a coordinate map.

I've been studying geodesics (not with much luck) and reading up on equations for the impact parameter and deflection angle in the schwarzschild metric which I have as :
$$\Delta\phi = \int^{r_{observed}}_{r_{omitted}} \stackrel{dr}{r\sqrt{r^{2}/b^{2} - 1 + R_{s}/r}}$$

Last edited: Oct 25, 2009
8. Oct 25, 2009

### A.T.

9. Oct 25, 2009

### nocks

Last edited by a moderator: May 4, 2017
10. Oct 25, 2009

### George Jones

Staff Emeritus
11. Oct 25, 2009

### nocks

Yes something like this is exactly what i'm aiming for at the moment. I noticed you focused on timelike geodesics for the app. Am I right in thinking that I would be using null geodesics for photon trajectories and timelike for the trajectory of an observer.

12. Oct 25, 2009

### George Jones

Staff Emeritus
Actually, my Exercise 6 is about lightlike (null) geodesics, which is what you want.

13. Oct 25, 2009

### Dmitry67

George, it is wonderful.
Do you know any animations about what observer would see inside the second horizon, near loop singularity?

14. Oct 26, 2009

### nocks

Would I still use geodesics for photon paths in polar coordinates?

15. Nov 1, 2009

### nocks

Could you give me some more information on plotting the trajectory on a coordinate map? I've been toying with equations for a while now and not making much progress
Thanks

16. Nov 1, 2009

### George Jones

Staff Emeritus
Yes, exactly. The geodesic equation when $\theta = \pi /2$ then is

$$\begin{equation*} \begin{split} \frac{d \phi}{d \lambda} &= \frac{L}{r^2} \\ \left( \frac{dr}{d \lambda} \right)^2 &= E^2 - L^2 W \left( r \left( \lambda \right) \right), \\ \end{split} \end{equation*}$$

where $W \left( r \right)$ is function that I'll specify later, and $E$ and $L$ are constants of motion.

In its present form, the second equation is a little difficult to implement on a computer since sometimes a positive square root is needed (increasing $r$) and sometimes a negative square root is needed (decreasing $r$) in the same photon orbit. To get aorund this, differentiate the second equation with respect to the affine parameter $\lambda$ taking into account that the $r$ in $W \left( r \right)$ is itself a function of $\lambda$.

What do you get for the second equation after this differentiation?

Last edited: Nov 21, 2009
17. Nov 1, 2009

### nocks

Given it's been a while since i've done any differentiation, is it simply:

$$\frac{dr}{d \lambda} \right) &= \frac{ E - L}{W \left( r \right)}$$

18. Nov 4, 2009

### nocks

what does W(r) define?

19. Nov 4, 2009

### Nabeshin

W(r) should be the effective potential of the system.

20. Nov 20, 2009

### nocks

Could anyone expand on this please? I would appreciate the help.
I have the the effective potential in the schwarzschild metric as (L being angular momentum)
$$V_{eff} = ( 1 - \frac{r_{s}}{r})(mc^{2} + \frac{L^{2}}{mr^{2}})$$

Would this be enough information to solve for r and $$\Phi$$ so that I could plot the trajectories.

Also could I use the same equation for plotting the course of an observer descending into the black hole?

21. Nov 20, 2009

### bcrowell

Staff Emeritus
Here is some python code I wrote that does pretty much what you're talking about. It's meant to be short and easy to understand, so it uses a pretty crude method of doing the numerical integration. If you want to do accurate numerical calculations of geodesics, you'd want to substitute a better integration method. There are various general-purpose subroutines out there, e.g., in the book Numerical Recipes in C. What my code does is to calculate the deflection of a light ray that grazes the sun. Actually, it calculates half he deflection for a ray that grazes the sun, with the mass of the sun scaled up by a factor of 1000 in order to keep the result from being overwhelmed by rounding errors in my el-cheapo integration method.

Code (Text):
import math

# constants, in SI units:
G = 6.67e-11         # gravitational constant
c = 3.00e8           # speed of light
m_kg = 1.99e30       # mass of sun
r_m = 6.96e8         # radius of sun

# From now on, all calculations are in units of the

# mass of sun, in units of the radius of the sun:
m_sun = (G/c**2)*(m_kg/r_m)
m = 1000.*m_sun

# Start at point of closest approach.
# initial position:
t=0
r=1 # closest approach, grazing the sun's surface
phi=-math.pi/2
# initial derivatives of coordinates w.r.t. lambda
vr = 0
vt = 1
vphi = math.sqrt((1.-2.*m/r)/r**2)*vt # gives ds=0, lightlike

l = 0    # affine parameter lambda
l_max = 20000.
epsilon = 1e-6 # controls how fast lambda varies
while l<l_max:
dl = epsilon*(1.+r**2) # giant steps when farther out
l = l+dl
# Christoffel symbols:
Gttr = m/(r**2-2*m*r)
Grtt = m/r**2-2*m**2/r**3
Grrr = -m/(r**2-2*m*r)
Grphiphi = -r+2*m
Gphirphi = 1/r
# second derivatives:
#  The factors of 2 are because we have, e.g., G^a_{bc}=G^a_{cb}
at   = -2.*Gttr*vt*vr
ar   = -(Grtt*vt*vt + Grrr*vr*vr + Grphiphi*vphi*vphi)
aphi = -2.*Gphirphi*vr*vphi
# update velocity:
vt = vt + dl*at
vr = vr + dl*ar
vphi = vphi + dl*aphi
# update position:
r = r + vr*dl
t = t + vt*dl
phi = phi + vphi*dl

# Direction of propagation, approximated in asymptotically flat coords.
# First, differentiate (x,y)=(r cos phi,r sin phi) to get vx and vy:
vx = vr*math.cos(phi)-r*math.sin(phi)*vphi
vy = vr*math.sin(phi)+r*math.cos(phi)*vphi
prop = math.atan2(vy,vx) # inverse tan of vy/vx, in the proper quadrant
prop_sec = prop*180.*3600/math.pi
print "final direction of propagation = %6.2f arc-seconds" % prop_sec

22. Nov 20, 2009

### nocks

I actually have the numerical recipes book next to me although I may avoid solving the elliptic integral https://www.physicsforums.com/showpost.php?p=2409536&postcount=7", and just use the approximation for light deflection, i.e. 4GM/bc$$^{2}$$, to get the einstein ring effect, and focus on the trajectory of the observer descending into the black hole.

Last edited by a moderator: Apr 24, 2017
23. Nov 20, 2009

### bcrowell

Staff Emeritus
Keep in mind that 4GM/bc2 is only a weak-field approximation. It won't give you the right answer if you're close to the black hole.

Last edited by a moderator: Apr 24, 2017
24. Nov 21, 2009

### George Jones

Staff Emeritus
No, this isn't quite right. Let's start with the left side. What is

$$\frac{d}{d\lambda} \left[ \left( \frac{dr}{d \lambda} \right)^2 \right] ?$$

Last edited: Nov 21, 2009
25. Nov 21, 2009

### nocks

Given

$$\left( \frac{dr}{d \lambda} \right)^2 &= E^2 - V^2(r) \right)$$
and
$$V^2(r) = \left(1 - \frac{2M}{r} \right)\frac{L^2}{r^2}$$

I have $$\frac{d^2r}{d\lambda^2} = -\frac{1}{2}\frac{d}{dr}V^2(r)$$