Help with using the Runge-Kutta 4th order method

  • Thread starter quasarLie
  • Start date
  • #1
51
0

Homework Statement


I have this equations of motion, I have this equations of motion for a schwarchild black hole, I wish to use the 4th order Runge-Kutta method to solve them for a body falling to the black hole from a distance r0 and with L = 0. My problem is I am struggling to apply this method to my system of ODE's so that I can program a method that can solve any system of ODE's using the formulas for RK4, I would like for someone to please run through one step of the method, so I can understand it better Thank

Homework Equations


[/B]
\begin{eqnarray}
\frac{dt}{d\tau} &=& (1- \frac{1}{r})^{-1} \frac{E}{m c^{2}}\\
\frac{d\phi}{d\tau} &= & \dfrac{L}{mr^{2}} \\
(\dfrac{dr}{d\tau})^{2} &=& \frac{E^{2}}{m^{2} c^{4}}(1- \frac{R_{s}}{r}) \dfrac{L^{2}}{m^{2} c^{2} r^{2}}
\end{eqnarray}​

The Attempt at a Solution


Now I know that for a general 1st order ODE's the 4th order Runge-Kutta formula's:
\begin{equation}
y_{i+1}=y_i + \frac{1}{6}(k_1+2k_2+2k_3+k_4)
\end{equation}
with
\begin{equation}
k_1 = f(t_n,y_n) \\ k_2 = f(t_n+\frac{1}{2}h,y_n+\frac{h}{2}k_1) \\ k_3 = f(t_n+\frac{1}{2}h,y_n+\frac{h}{2}k_2) \\ k_4 = f(t_n+\frac{1}{2}h,y_n+\frac{h}{2}k_3)
\end{equation}
 
Last edited by a moderator:

Answers and Replies

  • #2
34,158
5,778
You posted this originally in the Precalculus section -- this is hardly a precalc type of problem

Homework Statement


I have this equations of motion, I have this equations of motion for a schwarchild black hole, I wish to use the 4th order Runge-Kutta method to solve them for a body falling to the black hole from a distance r0 and with L = 0.
To the best of my knowledge, you can apply this method only to a single differential equation, not a system of differential equations.
quasarLie said:
My problem is I am struggling to apply this method to my system of ODE's so that I can program a method that can solve any system of ODE's using the formulas for RK4, I would like for someone to please run through one step of the method, so I can understand it better Thank
Which step are you uncertain of?
quasarLie said:

Homework Equations


[/B]
\begin{eqnarray}
\frac{dt}{d\tau} &=& (1- \frac{1}{r})^{-1} \frac{E}{m c^{2}}\\
\frac{d\phi}{d\tau} &= & \dfrac{L}{mr^{2}} \\
(\dfrac{dr}{d\tau})^{2} &=& \frac{E^{2}}{m^{2} c^{4}}(1- \frac{R_{s}}{r}) \dfrac{L^{2}}{m^{2} c^{2} r^{2}}
\end{eqnarray}​

The Attempt at a Solution


Now I know that for a general 1st order ODE's the 4th order Runge-Kutta formula's:
\begin{equation}
y_{i+1}=y_i + \frac{1}{6}(k_1+2k_2+2k_3+k_4)
\end{equation}
with
\begin{equation}
k_1 = f(t_n,y_n) \\ k_2 = f(t_n+\frac{1}{2}h,y_n+\frac{h}{2}k_1) \\ k_3 = f(t_n+\frac{1}{2}h,y_n+\frac{h}{2}k_2) \\ k_4 = f(t_n+\frac{1}{2}h,y_n+\frac{h}{2}k_3)
\end{equation}
 
  • #3
51
0
is this correct for the first equation
\begin{equation}
k_1=f(t_n,\tau_n) = (1- \frac{1}{r_n})^{-1} \frac{E}{m c^{2}}
\end{equation}
 
Last edited:
  • #4
phyzguy
Science Advisor
4,620
1,572
To the best of my knowledge, you can apply this method only to a single differential equation, not a system of differential equations.
This is not correct. RK4 works fine for systems of ODEs. You just have to interpret the [itex]y_n[/itex] as multicomponent vectors. In your case, [itex]\tau[/itex] is your independent variable, and your [itex] y_n[/itex] s are:
[tex] y_n = (t_n(\tau), \phi_n(\tau), r_n(\tau))[/tex]
Here is a link which explains in more detail.
 
  • Like
Likes scottdave
  • #5
34,158
5,778
To the best of my knowledge, you can apply this method only to a single differential equation, not a system of differential equations.
This is not correct. RK4 works fine for systems of ODEs. You just have to interpret the yn as multicomponent vectors.
I had something like that in mind, but didn't say it. What I was trying to get across was that you couldn't just plug the whole system of equations into RK. You're solving each equation for yi separately.
 
  • #6
51
0
I had something like that in mind, but didn't say it. What I was trying to get across was that you couldn't just plug the whole system of equations into RK. You're solving each equation for yi separately.
Yes i want to solve each esuation for y
 
  • #8
51
0
From post #2:
I want to know if this is correct for the second equation
\begin{equation}
k_{1}= \frac{L}{mr^{2}}\\
K_{2}= f(\phi(i)+\dfrac{h}{2},\dfrac{L}{mr_{i}^{2}}+ \frac{h}{2} \dfrac{L}{mr_{i}^{2}}\\
K_{3}= f(\phi_{i}+\frac{h}{2}, \dfrac{L}{mr_{i}^{2}}+\dfrac{h}{2}\dfrac{L}{mr_{i}^{2}})\\
K_{4}= f(\phi_{i}+h,\dfrac{L}{mr_{i}^{2}}+h\dfrac{L}{mr_{i}^{2}})
\end{equation}
 
  • #9
34,158
5,778
From post #1:
## \frac{dt}{d\tau} = (1- \frac{1}{r})^{-1} \frac{E}{m c^{2}} ##
## \frac{d\phi}{d\tau} = \frac{L}{mr^2} ##
## (\frac{dr}{d\tau})^{2} = \frac{E^{2}}{m^{2} c^{4}}(1- \frac{R_{s}}{r}) \frac{L^{2}}{m^{2} c^{2} r^{2}}##
I think you're going to have to do some work on the third equation before you can tackle the other two equations, since neither of the first two equations actually involves τ.
The third equation could be rewritten as ##\frac {dr}{d\tau} = \frac E {mc^2} \sqrt{1 - \frac {R_s}r} \frac L {mcr}##. Here I'm assuming that all of the squared quantities are nonnegative, so for instance, ##\sqrt{m^2} = m##.
If you solve this (manually), you will have r as a function of ##\tau##, and you can replace r in the first and second equations.
Then the first equation will have dt/dτ equal to an expression that actually involves τ, and similar for the second equation.
 
  • #10
51
0
From post #1:
## \frac{dt}{d\tau} = (1- \frac{1}{r})^{-1} \frac{E}{m c^{2}} ##
## \frac{d\phi}{d\tau} = \frac{L}{mr^2} ##
## (\frac{dr}{d\tau})^{2} = \frac{E^{2}}{m^{2} c^{4}}(1- \frac{R_{s}}{r}) \frac{L^{2}}{m^{2} c^{2} r^{2}}##
I think you're going to have to do some work on the third equation before you can tackle the other two equations, since neither of the first two equations actually involves τ.
The third equation could be rewritten as ##\frac {dr}{d\tau} = \frac E {mc^2} \sqrt{1 - \frac {R_s}r} \frac L {mcr}##. Here I'm assuming that all of the squared quantities are nonnegative, so for instance, ##\sqrt{m^2} = m##.
If you solve this (manually), you will have r as a function of ##\tau##, and you can replace r in the first and second equations.
Then the first equation will have dt/dτ equal to an expression that actually involves τ, and similar for the second equation.
Yes I have work to do with last equation to apply the RK4, i just want ti know if what i did is correct to use it in a C++ program to solve the quations. Thank you for your answers
 
  • #11
34,158
5,778
Yes I have work to do with last equation to apply the RK4, i just want ti know if what i did is correct to use it in a C++ program to solve the quations. Thank you for your answers
No, it doesn't look correct to me.
Here is your second equation: ## \frac{d\phi}{d\tau} = \frac{L}{mr^2} ##
To use Runge-Kutta, the equation needs to be in this form: ##\frac {d\phi}{d\tau} = f(\phi, \tau)##. That is, the right side of the equation has to be in terms of ##\phi## and ##\tau##. What you have in your second equation immediately above doesn't involve either of these variables. That's why I'm suggesting that you find a solution to the third equation, so that you can eventually get r in terms of ##\tau##, then substitute for r in the first two equations.
 
Last edited:
  • #12
51
0
No, it doesn't look correct to me.
Here is your second equation: ## \frac{d\phi}{d\tau} = \frac{L}{mr^2} ##
To use Runge-Kutta, the equation needs to be in this form: ##\frac {d\phi}{d\tau} = f(\phi, \tau)##. That is, the right side of the equation has to be in terms of ##\phi## and ##\tau##. What you have in your second equation immediately above doesn't involve either of these variables. That's why I'm suggesting that you find a solution to the third equation, so that you can eventually get r in terms of ##\tau##, then substitute for r in the first two equations.
This is what I did for the RK2. Do you think that it's false ? (it's a python program)
tRK2[i+1] = t + dtau * (t + dtau/2 + (1/(1-(Rs/R)) * E) + (dtau/2)*(1/(1-(Rs/R)) * E))
phiRK2[i+1] = phi + dtau * (phi + dtau/2 + (L/(m*R**2)) + (dtau/2)*(L/(m*R**2)))
RRK2[i+1] = R + dtau * (R + dtau/2 + V + (dtau/2)*V)
VRK2[i+1] = VRK2 - dtau * (VRK2 + dtau/2 + ((Rs/(2*R**2)) + (3/2)*((Rs*L**2 - R*L**2)/(R**4))) + (dtau/2)*((Rs/(2*R**2)) + (3/2)*((Rs*L**2 - R*L**2)/(R**4))))
 
  • #13
34,158
5,778
This is what I did for the RK2. Do you think that it's false ? (it's a python program)
Python:
 tRK2[i+1] = t[I] + dtau * (t[I] + dtau/2 + (1/(1-(Rs/R[I])) * E) + (dtau/2)*(1/(1-(Rs/R[I])) * E))
        phiRK2[i+1] = phi[I] + dtau * (phi[I] + dtau/2 + (L/(m*R[I]**2)) + (dtau/2)*(L/(m*R[I]**2)))
        RRK2[i+1] = R[I] + dtau * (R[I] + dtau/2 + V[I] + (dtau/2)*V[I])
        VRK2[i+1] = VRK2[I] - dtau * (VRK2[I] + dtau/2 + ((Rs/(2*R[I]**2)) + (3/2)*((Rs*L**2 - R[I]*L**2)/(R[I]**4))) + (dtau/2)*((Rs/(2*R[I]**2)) + (3/2)*((Rs*L**2 - R[I]*L**2)/(R[i]**4))))[/i][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/I][/QUOTE]
Your program is unreadable -- all expressions such as t[i] get converted to BBCode for italics and disappear. Also, when you post code, you should use code tags -- see https://www.physicsforums.com/threa...-for-nearly-200-programming-languages.774303/ for more information. Using code tags retains the indentation, which is critical for Python code.
Second, I've invested a fair amount of time in perusing what you've done for RK4. Why don't you proceed with that, rather than persue a different direction?
 
  • #14
51
0
Your program is unreadable -- all expressions such as t[i] get converted to BBCode for italics and disappear. Also, when you post code, you should use code tags -- see https://www.physicsforums.com/threa...-for-nearly-200-programming-languages.774303/ for more information. Using code tags retains the indentation, which is critical for Python code.
Second, I've invested a fair amount of time in perusing what you've done for RK4. Why don't you proceed with that, rather than persue a different direction?
Thank you for the information. and I just want to show you what i did for RK2 (which is correct) because I want to do the same thing with RK4 this is why I wrote the K_i
 
  • #15
Ray Vickson
Science Advisor
Homework Helper
Dearly Missed
10,706
1,728
Yes I have work to do with last equation to apply the RK4, i just want ti know if what i did is correct to use it in a C++ program to solve the quations. Thank you for your answers
Besides what Mark44 said in post #9, another possibiity is that
$$ \frac {dr}{d\tau} = -\frac E {mc^2} \sqrt{1 - \frac {R_s}r} \frac L {mcr}$$
That will give different dynamics from the choice
$$ \frac {dr}{d\tau} = \frac E {mc^2} \sqrt{1 - \frac {R_s}r} \frac L {mcr},$$
so you need to try to figure out which is the relevant form.
 

Related Threads on Help with using the Runge-Kutta 4th order method

  • Last Post
Replies
2
Views
2K
Replies
0
Views
8K
Replies
6
Views
1K
Replies
2
Views
2K
  • Last Post
Replies
2
Views
2K
Replies
1
Views
3K
  • Last Post
Replies
1
Views
4K
Replies
5
Views
2K
  • Last Post
Replies
8
Views
2K
Top