Using ode45 to solve a system of second order differential equations.

Hello, I am trying to solve an orbit problem using the J2 disturbance. If you do not know what that is, it is irrelevant anyways. All I need to know is how to numerically integrate a system of second order differential equations.

An example of one of them is:

x'' = -mu*x/(x^2+y^2+z^2)^(3/2)*(1 - 3/2*J2*(Re/(x^2+y^2+z^2)^(1/2))^2*(5*z^2/(x^2+y^2+z^2)-1))

mu, Re, and J2 are all constants.

This and the other two equations, y'' and z'', which are similar, define the radius of a satellite while the first derivatives define the velocity in Cartesian coordinates.

When initializing my function into ode45, I tried to use a matrix, not a one dimensional array, and it would not let me. I am not sure if I can use symbolic variables to do a nested ode45 command to solve say, z'' in terms of x and y, and then y'' and then x'' to find the solution for that timestep and repeat.

If anyone has any hints, that would be very useful. Thanks!

Answer 1

[Torsten]

The second and third example under

http://de.mathworks.com/help/matlab/ref/ode45.html

are second-order equations solved with ODE45.

In your case, solve all equations simultaneously (set y(1)=x, y(2)=x', y(3)=y, y(4)=y', y(5)=z, y(6)=z').

Use a one-dimensional array for initialization of x, x', y, y', z and z'.

Don't use symbolic algebra.

Best wishes

Torsten.

[Robert]

I used that resource to the best that I could; however, the only system of equations example it uses is first order. I may be confused by this, but if you use what you said, wouldn't the code be assuming that your y's go up to the 5th derivative of the variable in question, not zeroth and first three times?

I tried that approach but it did not give me the answer that I know to be correct. There is a chance that I could have inputted something incorrectly, but is my theory about the 5th derivative incorrect?

[Torsten]

Input the system

y1'=y2
y2'=-mu*y1/(y1^2+y3^2+y5^2)^(3/2)*(1 - 3/2*J2*(Re/(y1^2+y3^2+y5^2)^(1/2))^2*(5*y5^2/(y1^2+y3^2+y5^2)-1))
y3'=y4
y4'=expression for y''
y5'=y6
y6'=expression for z''

and you'll get x as y1, x' as y2, y as y3, y' as y4, z as y(5) and z' as y6.

If it does not give the solution you expect, please include your code.

Best wishes

Torsten.