Solving Orbital Motion using ODE45 - MY Problem

[aerospace37]

MY problom is easy,but i can't get the reasonable answer.problom dipicted:
In the Cartesian coordinates,A satellite'higth is 200KM,and just in two dimensions.the initial position is （-6571KM，0），and initial velocity is （0，-7.8KM/S）.I need the final running solution of program.I want to solve the problom in ode45.The result figure just is a circle
program as following:
============================================
function Yd=orbit(t,y)
global u
rx=y(1);
ry=y(2);
vx=y(3);
vy=y(4);
rr=sqrt(rx^2+ry^2);
vxy=[vx;vy];
xy=[rx;ry];
Yd=[vxy;-u*xy/rr^3];
===================================

function wlow3
global u
u=3.986e14;
t0=0;tf=24*60*60;
tspan=[t0,tf]; %
y0=[-6.571e6;0;0;-7.8e3];%
[t,YY]=ode45('orbit',tspan,y0);
X=YY(:,1);
Y=YY(:,2);
plot(X,Y);
xlabel('x')
ylabel('y')
hold on
axis('image')

==========================================================
I hope anyone who researches orbit or others can communicates with me.
thans a lot.

 

[jvicens]

Thank you for sharing your problem with us. I am a scientist who specializes in orbital dynamics and I would be happy to offer some help with your program.

Firstly, it is important to note that the final running solution of your program will depend on the specific values of the initial conditions and the gravitational parameter (u) that you have chosen. In your program, the gravitational parameter is set to u=3.986e14, which is the standard gravitational parameter for Earth. However, the initial conditions you have provided do not seem to correspond to a circular orbit around Earth.

To solve this problem using ode45, you will need to define the equations of motion for a satellite in Cartesian coordinates. These equations can be derived from Newton's laws of motion and the law of gravitation. In your program, the orbit function is defined as Yd=[vxy;-u*xy/rr^3], which is the correct form for the equations of motion in Cartesian coordinates. However, there are a few changes that need to be made in order to ensure that the program runs correctly.

Firstly, you have defined the initial conditions as y0=[-6.571e6;0;0;-7.8e3], which corresponds to an initial position of (-6571 km, 0) and an initial velocity of (0, -7.8 km/s). However, the units for position and velocity should be consistent. In this case, the units for position are in meters and the units for velocity are in meters per second. Therefore, the initial conditions should be defined as y0=[-6571000;0;0;-7800].

Secondly, the equations of motion should be written in terms of the state variables x and y, rather than rx and ry. This can be done by modifying the orbit function to the following:

function Yd=orbit(t,y)
global u
x=y(1);
y=y(2);
vx=y(3);
vy=y(4);
rr=sqrt(x^2+y^2);
vxy=[vx;vy];
xy=[x;y];
Yd=[vxy;-u*xy/rr^3];

Finally, to plot the orbit, you will need to convert the state variables x and y back to kilometers, as follows:

X=YY(:,1)/1000;
Y=YY(:,2)/1000;
plot(X,Y);
xlabel('x (km)')