Couple projectile motion equations function
2 ビュー (過去 30 日間)
古いコメントを表示
I am creating a function to solve these coupled equations, where all constants are known:
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/288418/image.png)
I am currently using the code below, but I am not getting the correct results. All of the constants are specified in my code, b If anyone could offer any help that would be great!
% Set constants for different projectiles (baseball, football, soccer ball)
global m cd A rho g theta v0 x0 y0
m= [0.145 0.42 0.43]; % mass of projectile [kg]
cd= [0.3 0.5 0.25]; % drag coefficients
A= [0.00426 0.02318 0.038]; % cross sectional area of projectile [m^2]
V= [0.0002194 0.00462115 0.00573547]; % volume of projectile
rho= [m(1)/V(1) m(2)/V(2) m(3)/V(3)]; % density
g= 9.8; % gravity [m/s^2]
theta= 45; % initial launch angle [rad]
x0= [0 0 0]; % initial horizontal position
y0= [1.8288 1.8288 0]; % initial vertical position
v0= [18 18 22]; % initial velocity (m/s)
%set projectile, 1=baseball, 2=football, 3=soccer ball
proj=2;
m= m(proj);
cd= cd(proj);
A= A(proj);
rho= rho(proj);
x0= x0(proj);
y0= y0(proj);
v0= v0(proj);
% call function
t_end= 2*v0*sin(theta)/g; % time until projectile reaches the ground
[t y]= ode45(@my_fun,[0 t_end],[x0 y0 v0 theta]);
% plotting results
figure(1)
plot(t,y(:,1))
grid on
xlabel('time [s]')
ylabel('horizontal position [m]')
title('plot of horizontal location vs time')
figure(2)
plot(t,y(:,2))
grid on
xlabel('time [s]')
ylabel('vertical position [m]')
title('plot of vertical location vs time')
figure(3)
plot(y(:,1),y(:,2))
grid on
xlabel('horizontal position [m]')
ylabel('vertical position [m]')
title('plot of vertical location vs horizontal position')
%% Function
function dy=my_fun(t,y) % y is a matrix containing [x y v theta]
global m cd A rho g theta x0 y0 v0
dy=zeros(4,1);
dy(1)= y(3); %dx/dt
dy(2)= y(3); %dy/dt
dy(3)= -(cd*A*rho)/(2*m)*sqrt((dy(1))^2+(dy(2))^2)*dy(1); % d2x/dt2
dy(4)= -g-(cd*A*rho)/(2*m)*sqrt((dy(1))^2+(dy(2))^2)*dy(2); % d2y/dt2
end
0 件のコメント
採用された回答
Ameer Hamza
2020 年 4 月 28 日
In your my_fun you have a line
dy(2)= y(3); %dy/dt
It should be
dy(2)= y(4); %dy/dt
Everything else is correct according to the equation in your question.
2 件のコメント
Nicholas Moose
2020 年 12 月 14 日
It looks like y(4) is a launch angle. Shouldn't it be something like dy(1) = y(3)*cos(y(4)) and dy(2) = y(3)*sin(y(4))?
James Tursa
2020 年 12 月 14 日
編集済み: James Tursa
2020 年 12 月 14 日
No. There are four states: x, y, vx, vy. These are defined to be y(1), y(2), y(3), and y(4) in the code. The dx/dt value is always going to be vx, which is simply y(3), and the dy/dt value is always going to be vy, which is simply y(4).
The first dx/dt and dy/dt values may in fact correspond to a launch angle, but that is only for the initial point.
That being said, the initial conditions that are passed to ode45( ) do not involve theta directly as was written. They should only be [x0 y0 vx0 vy0] where all of these values are scalars. So in this case since we are working with launch angle theta it would be:
[x0 y0 v0*cosd(theta) v0*sind(theta)]
Also this
theta= 45; % initial launch angle [rad]
should be this
theta= 45; % initial launch angle [deg]
and downstream in the code you should be using sind(theta) instead of sin(theta).
その他の回答 (0 件)
参考
製品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!