How to use ODE 45 to integrate equations of motion?

10 vues (au cours des 30 derniers jours)
Stephen Mixon
Stephen Mixon le 28 Oct 2019
Commenté : Cyrus Tirband le 29 Oct 2019
I am modeling a 3d projectile and need help using the ode45 command. I need it to integrate the equations of motion based on my simulation and for it to end when z=0. I need to take this numerical results and sumperimpose them in a plot with my previous data. I need to plot the errors in my data beween my parameter set and this ode 45 set.
This is my code for the simulation
figure('name','Projectile Motion')
xlabel('X (m)')
ylabel('Y (m)')
zlabel('Z (m)')
% This sets the viewing angle
view([1 -1 1])
% This sets the limits on the x- and y-axes
xlim([0,3000])
ylim([0,3000])
zlim([0,1000])
box('on')
grid('on')
hold on
% az is azimuth angle, el is elevation
V0 = 150;
az = 45;
el= 30;
t = [0:0.01:25]';
vx0 = V0*cosd(az)*cosd(el);
vy0 = V0*sind(az)*cosd(el);
vz0 = V0*sind(el);
x = vx0.*t;
y = vy0.*t;
z = vz0*t-0.5*9.81.*t.^2;
[~,ind]=min(abs(z));
plot3(x,y,z)
%Time of flight
T =(2*vz0)/9.81
  7 commentaires
Stephen Mixon
Stephen Mixon le 28 Oct 2019
Write MATLAB code that uses the ode45 package to numerically integrate the equations of
motion and simulate the flight of the projectile. The integration should stop when the projectile
reaches . The simulation and subsequent post-processing should give you all the
z=0
parameters that you obtained analytically in Part 1a.
Compare the numerical results to the analytical results for each parameter computed in Part 1a
by superimposing plots. Time should be on the x-axis and the parameter on the y-axis. Use
dashed lines for analytical results and a continuous line of a different color for the numerical
results.
darova
darova le 28 Oct 2019
What should be the result of integration? What is the parameter? t?

Connectez-vous pour commenter.

Réponses (1)

Cyrus Tirband
Cyrus Tirband le 28 Oct 2019
First, you need to write down your equations of motion. In your simple case, you have:
with initial conditions
Of course, these integrals are fairly easy to solve and you don't really need to verify these with a numerical solver, but for your assignment you are asked to compare them.
Knowing the equations, the ode function becomes easy to set up. Each of them are independent, and almost the same, so
function dzdt = odefunc(t,z,g)
dzdt = [y(2), -g];
end
similar for x, y (or set g to zero)
then call the solver
[t,z] = ode45(@(t,z), odefunc(t,z,g), [0 tend],[0; vz0]);
The solutions can be seen to overlap perfectly.
ode_example.png
  2 commentaires
Stephen Mixon
Stephen Mixon le 28 Oct 2019
I used the code you gave me and inseted inital conditions to 0 and got this error
Function definitions in a script must appear at the end of the file.
Move all statements after the "odefunc" function definition to before the first local
function definition.
here is the updated code
figure('name','Projectile Motion')
xlabel('X (m)')
ylabel('Y (m)')
zlabel('Z (m)')
% This sets the viewing angle
view([1 -1 1])
% This sets the limits on the x- and y-axes
xlim([0,3000])
ylim([0,3000])
zlim([0,1000])
box('on')
grid('on')
hold on
% az is azimuth angle, el is elevation
V0 = 150;
az = 45;
el= 30;
t = [0:0.01:25]';
vx0 = V0*cosd(az)*cosd(el);
vy0 = V0*sind(az)*cosd(el);
vz0 = V0*sind(el);
x = vx0.*t;
y = vy0.*t;
z = vz0*t-0.5*9.81.*t.^2;
[~,ind]=min(abs(z));
plot3(x,y,z)
%Time of flight
T =(2*vz0)/9.81;
x_a=0;
y_a=0;
z_a=0;
x_v=vx0;
x_v=vy0;
z_v=vz0;
x0=0;
y0=0;
z0=0;
function dzdt = odefunc(t,z,g)
dzdt = [y(2), -g];
end
[t,z] = ode45(@(t,z), odefunc(t,z,g), [0 tend],[0; vz0]);
Cyrus Tirband
Cyrus Tirband le 29 Oct 2019
The error message is fairly unambiguous. You cannot define functions in the middle of the script, you have to put them at the end.
If you're using MATLAB 2016a or earlier, you'll have to define them in their own file, or make the main script a function.
Additionally, you'll have to define a few parameters before my snippets will generate the solution.

Connectez-vous pour commenter.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by