Projectile motion of a cannonball with varying drag
11 次查看(过去 30 天)
显示 更早的评论
Hey everyone,
https://imgur.com/a/yUTmYXh
For a project, I need to simulate the projectile motion of a "Paris Cannon". I need to consider the air resistance changing according to the altitude with the following equation:
where y0 = 1000m.
And the drag is:
So far, my code looks like this:
clc
clear
m=94;
A=pi*(.21)^2;
C=.12;
p0= 1.225;
g=9.81;
%Initial Conditions
dt= .001;
x(1)=0;
y(1)=0;
v=10000;
angle=20 ;
vx=v*cosd(angle);
vy=v*sind(angle);
t(1)=0 ;
i=1;
y0 = 1000;
while min(y)> -.001
p = p0*exp(-y/y0);
Drag = p*C*A/2*v^2;
ax=-(Drag/m)*(vx*vx+vy*vy)^0.5*vx;
ay=-g-(Drag/m)*(vx*vx+vy*vy)^0.5*vy;
vx=vx+ax*dt;
vy=vy+ay*dt;
x(i+1)=x(i)+vx*dt+.5*ax*dt^2;
y(i+1)=y(i)+vy*dt+.5*ay*dt^2;
t(i+1)=t(i)+dt;
i=i+1;
end
plot(x,y,'-')
xlabel('x distance (m)')
ylabel('y distance (m)')
pause(0.001);
I cannot get the correct simulation here. If I take drag as constant and place it in initial conditions, it works fine. But when I put it in while loop, the whole simulation gives me absolute nonsense. I would appriciate your help.
Thanks.
0 个评论
回答(1 个)
Michael Madelaire
2018-12-30
Hi,
I am not sure about the equations you use. But here is my attempt.
Try running it, and if it is what you are looking for, maybe I can explain how I did it.
m=94;
A=pi*(.21)^2;
C=.12;
p0= 1.225;
g=9.81;
%Initial Conditions
dt= .001;
x=0;
y=0;
v=10000;
angle=20 ;
vx=v*cosd(angle);
vy=v*sind(angle);
t=0 ;
y0 = 1000;
i=1;
% Half step
p = p0*exp(-y/y0);
D = C*(p*v^2)/2*A;
ax = -D/m;
ay = -g;
while y(i) >= 0
i = i+1;
vx(i) = vx(i-1)+ax(i-1)*dt;
vy(i) = vy(i-1)+ay(i-1)*dt;
x(i) = x(i-1)+(vx(i-1)+vx(i))/2*dt;
y(i) = y(i-1)+(vy(i-1)+vy(i))/2*dt;
p = p0*exp(-y(i)/y0);
v = sqrt(vx(i)^2+vy(i)^2);
D = C*(p*v^2)/2*A;
ax(i) = -D/m;
ay(i) = -g;
end
figure;
plot(x,y,'-')
xlabel('x distance (m)')
ylabel('y distance (m)')
5 个评论
Michael Madelaire
2018-12-30
Yes, but you need to find what differential equations you are solving and to what order of Runge-Kutta you want to use.
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!