light ray paths near schwarzschild blackhole

[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

[nocks]

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

[ZikZak]

You might start with this paper. (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 this paper. (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 'light deflection near neutron stars' (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

http://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

http://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 I mentioned earlier (http://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 I mentioned earlier (http://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/bc2 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