Light ray paths near schwarzschild blackhole

In summary, the conversation is about simulating light ray paths near black holes, specifically starting with Schwarzschild black holes and working towards Kerr-Newman black holes. The speaker is a computer science expert who is looking for advice on how to plot the trajectory of a light ray near a black hole on a coordinate map. They have been studying geodesics and equations for impact parameter and deflection angle, but have not made much progress. They have received helpful resources and advice from others in the conversation and are interested in plotting both photon and observer trajectories. The equations for the photon trajectory use the geodesic equation and the effective potential, while the equations for the observer trajectory use the effective potential and angular momentum.
  • #1
nocks
24
0
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
 
Last edited:
Physics news on Phys.org
  • #2
Anyone have any advice?
I'm currently going over lagrangian mechanics
 
  • #5
Last edited by a moderator:
  • #6
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
George Jones said:
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 :
[tex]\Delta\phi = \int^{r_{observed}}_{r_{omitted}} \stackrel{dr}{r\sqrt{r^{2}/b^{2} - 1 + R_{s}/r}} [/tex]
 
Last edited:
  • #9
Last edited by a moderator:
  • #11
George Jones said:
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.
 
  • #12
nocks said:
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.
 
  • #13
George, it is wonderful.
Do you know any animations about what observer would see inside the second horizon, near loop singularity?
 
  • #14
Would I still use geodesics for photon paths in polar coordinates?
 
  • #15
nocks said:
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
 
  • #16
nocks said:
Would I still use geodesics for photon paths in polar coordinates?

Yes, exactly. The geodesic equation when [itex]\theta = \pi /2[/itex] then is

[tex]\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*}
[/tex]

where [itex]W \left( r \right)[/itex] is function that I'll specify later, and [itex]E[/itex] and [itex]L[/itex] 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 [itex]r[/itex]) and sometimes a negative square root is needed (decreasing [itex]r[/itex]) in the same photon orbit. To get aorund this, differentiate the second equation with respect to the affine parameter [itex]\lambda[/itex] taking into account that the [itex]r[/itex] in [itex]W \left( r \right)[/itex] is itself a function of [itex]\lambda[/itex].

What do you get for the second equation after this differentiation?
 
Last edited:
  • #17
Given it's been a while since I've done any differentiation, is it simply:

[tex]\frac{dr}{d \lambda} \right) &= \frac{ E - L}{W \left( r \right)}[/tex]
 
  • #18
what does W(r) define?
 
  • #19
W(r) should be the effective potential of the system.
 
  • #20
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)
[tex]V_{eff} = ( 1 - \frac{r_{s}}{r})(mc^{2} + \frac{L^{2}}{mr^{2}})[/tex]

Would this be enough information to solve for r and [tex]\Phi[/tex] 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
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:
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
 
  • #22
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[tex]^{2}[/tex], to get the einstein ring effect, and focus on the trajectory of the observer descending into the black hole.
 
Last edited by a moderator:
  • #23
nocks said:
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[tex]^{2}[/tex], 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.
 
Last edited by a moderator:
  • #24
George Jones said:
Yes, exactly. The geodesic equation when [itex]\theta = \pi /2[/itex] then is

[tex]\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*}
[/tex]

where [itex]W \left( r \right)[/itex] is function that I'll specify later, and [itex]E[/itex] and [itex]L[/itex] 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 [itex]r[/itex]) and sometimes a negative square root is needed (decreasing [itex]r[/itex]) in the same photon orbit. To get aorund this, differentiate the second equation with respect to the affine parameter [itex]\lambda[/itex] taking into account that the [itex]r[/itex] in [itex]W \left( r \right)[/itex] is itself a function of [itex]\lambda[/itex].

What do you get for the second equation after this differentiation?
nocks said:
Given it's been a while since I've done any differentiation, is it simply:

[tex]\frac{dr}{d \lambda} \right) &= \frac{ E - L}{W \left( r \right)}[/tex]

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

[tex]\frac{d}{d\lambda} \left[ \left( \frac{dr}{d \lambda} \right)^2 \right] ?[/tex]
 
Last edited:
  • #25
George Jones said:
...
Given

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

I have [tex] \frac{d^2r}{d\lambda^2} = -\frac{1}{2}\frac{d}{dr}V^2(r)[/tex]
 
  • #26
nocks said:
[tex]V^2(r) = \left(1 - \frac{2M}{r} \right)\frac{L^2}{r^2}[/tex]
Guess I should state that

[tex]\frac{d}{dr}V^2(r) = \frac{2L^2M}{r^4} - \frac{2L^2(\left 1 - \frac{2M}{r} \right)}{r^3}[/tex]
 
  • #27
After simplifying the right side of
nocks said:
[tex]\frac{d}{dr}V^2(r) = \frac{2L^2M}{r^4} - \frac{2L^2(\left 1 - \frac{2M}{r} \right)}{r^3}[/tex]

and plugging back into
nocks said:
[tex] \frac{d^2r}{d\lambda^2} = -\frac{1}{2}\frac{d}{dr}V^2(r)[/tex]

what do you get?
 
Last edited:
  • #28
That would be

[tex]\frac{d^2r}{d\lambda^2} = - \frac{L^2(3M-r)}{r^4}[/tex]

Which gives me the closest radius for a stable orbit of a photon as 3M
 
  • #29
nocks said:
That would be

[tex]\frac{d^2r}{d\lambda^2} = - \frac{L^2(3M-r)}{r^4}[/tex]

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?
 
  • #30
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 [itex]\left( r , \phi \right)[/itex] position from which the photon is fired and the direction in which the photon is fired. These initial conditions determine [itex]L[/itex]. Alternatively, [itex]L[/itex] could be used as one of the initial conditions.

Much more on all of this later.
 
  • #31
nocks said:
That would be

[tex]\frac{d^2r}{d\lambda^2} = - \frac{L^2(3M-r)}{r^4}[/tex]

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 [itex]p = dr/d\lambda[/itex], so that [itex]dp/dt = d^2r/d\lambda^2[/itex]. The set of equations that describes the worldline of a photon then is

[tex]
\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*}
[/tex]

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.
 
  • #32
Having a bit of trouble with this but I have p as

p = [tex]\frac{L^2(2M-r)}{2r^3}[/tex]
 
  • #33
wow it's been a while since I looked at this project.
Thought i'd come back to it :)

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

[tex]
\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*}
[/tex]

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
 
Last edited:

Similar threads

  • Special and General Relativity
Replies
4
Views
288
  • Special and General Relativity
Replies
20
Views
686
  • Special and General Relativity
Replies
4
Views
1K
  • Special and General Relativity
Replies
8
Views
2K
  • Special and General Relativity
Replies
17
Views
2K
  • Special and General Relativity
Replies
5
Views
1K
Replies
58
Views
4K
  • MATLAB, Maple, Mathematica, LaTeX
Replies
3
Views
969
  • Special and General Relativity
Replies
2
Views
2K
  • Special and General Relativity
Replies
10
Views
1K
Back
Top