Light ray paths near schwarzschild blackhole

[nocks]

hey there, I'm interested in (eventually) simulating light ray paths near black holes, starting with schwarzschild black holes and working my way to kerr-newman black holes.
I have a good understanding of the nature of black holes 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

 

[nocks]

Anyone have any advice?
I'm currently going over lagrangian mechanics

 

[ZikZak]

You might start with http://lanl.arxiv.org/abs/astro-ph/9312003"

 

[Naty1]

Can't help much but this animation of gravitational lensing might be of interest...

http://en.wikipedia.org/wiki/Black_hole#History

 

[nocks]

You might start with http://lanl.arxiv.org/abs/astro-ph/9312003"

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

 

[George Jones]

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.

 

[nocks]

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.

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}}$$

 

[A.T.]

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.

Here some material.
http://www.spacetimetravel.org/reiseziel/reiseziel1.html

 

[nocks]

Here some material.
http://www.spacetimetravel.org/reiseziel/reiseziel1.html

I stumbled upon this last week, the paper on http://www.spacetimetravel.org/licht/light.pdf" by ute kraus has been quite usefull. Appreciate the help though :)

 

[George Jones]

for now I would like to simply plot the trajectory of several photons on a coordinate map.

You mean like Exercise 6 in the (imperfect) animation attached to

https://www.physicsforums.com/showthread.php?p=1091901#post1091901?

I won't have access to my notes and books until tomorrow.

 

[nocks]

You mean like Exercise 6 in the (imperfect) animation attached to

https://www.physicsforums.com/showthread.php?p=1091901#post1091901?

I won't have access to my notes and books until tomorrow.

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.

 

[George Jones]

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.

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

 

[Dmitry67]

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

 

[nocks]

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

 

[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.

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

 

[George Jones]

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

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?

 

[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)}$$

 

[nocks]

what does W(r) define?

 

[Nabeshin]

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

 

[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?

 

[bcrowell]

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.

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
# radius of the sun.

# 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

 

[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.

 

[bcrowell]

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.

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

 

[George Jones]

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?

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)}$$

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] ?$$

 

[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)$$

 

[nocks]

$$V^2(r) = \left(1 - \frac{2M}{r} \right)\frac{L^2}{r^2}$$

Guess I should state that

$$\frac{d}{dr}V^2(r) = \frac{2L^2M}{r^4} - \frac{2L^2(\left 1 - \frac{2M}{r} \right)}{r^3}$$

 

[George Jones]

After simplifying the right side of

$$\frac{d}{dr}V^2(r) = \frac{2L^2M}{r^4} - \frac{2L^2(\left 1 - \frac{2M}{r} \right)}{r^3}$$

and plugging back into

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

what do you get?

 

[nocks]

That would be

$$\frac{d^2r}{d\lambda^2} = - \frac{L^2(3M-r)}{r^4}$$

Which gives me the closest radius for a stable orbit of a photon as 3M

 

[nocks]

That would be

$$\frac{d^2r}{d\lambda^2} = - \frac{L^2(3M-r)}{r^4}$$

Which gives me the closest radius for a stable orbit of a photon as 3M

Apologies for the lack of basic physics understanding but for a photon on approach what would it's angular momentum be?

 

[George Jones]

We now have a system of three first-order coupled differential equations that determine the worldlines for "photons". In order to solve these equations, we need three initial conditions. Imagine firing a photon from a laser. One possibility for the three initial conditions: the $\left( r , \phi \right)$ position from which the photon is fired and the direction in which the photon is fired. These initial conditions determine $L$. Alternatively, $L$ could be used as one of the initial conditions.

Much more on all of this later.

 

[George Jones]

That would be

$$\frac{d^2r}{d\lambda^2} = - \frac{L^2(3M-r)}{r^4}$$

Which gives me the closest radius for a stable orbit of a photon as 3M

Reduce this second-order equation to two first-order equations by setting $p = dr/d\lambda$, so that $dp/dt = d^2r/d\lambda^2$. The set of equations that describes the worldline of a photon then is

$$\begin{equation*} \begin{split} \frac{d \phi}{d \lambda} &= \frac{L}{r^2} \\ \frac{dr}{d\lambda} &= p \\ \frac{dp}{d \lambda} &= \frac{L^2(r - 3M)}{r^4}.\\ \end{split} \end{equation*}$$

Assuming that all required values are given, write a few lines of (pseudo)code that uses the simplest, most intuitive method (Euler's method) to solve these equations.

 

[nocks]

Having a bit of trouble with this but I have p as

p = $$\frac{L^2(2M-r)}{2r^3}$$

 

[nocks]

wow it's been a while since I looked at this project.
Thought i'd come back to it :)

Reduce this second-order equation to two first-order equations by setting $p = dr/d\lambda$, so that $dp/dt = d^2r/d\lambda^2$. The set of equations that describes the worldline of a photon then is

$$\begin{equation*} \begin{split} \frac{d \phi}{d \lambda} &= \frac{L}{r^2} \\ \frac{dr}{d\lambda} &= p \\ \frac{dp}{d \lambda} &= \frac{L^2(r - 3M)}{r^4}.\\ \end{split} \end{equation*}$$

Assuming that all required values are given, write a few lines of (pseudo)code that uses the simplest, most intuitive method (Euler's method) to solve these equations.

A few questions. Since I lack a good maths background, I would appreciate some advice on how I would implement these 3 equations into a numerical solver like Euler's or Runge-Kutta to get back some results that I could plot.
e.g. Do I pre-define L?, how do I know which direction the photon is travelling?

currently trying to use this : http://www.ee.ucl.ac.uk/~mflanaga/java/RungeKutta.html