Projectile motion of a cannonball with varying drag

Question

Ozan Terzioglu 2018-12-30

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/437689-projectile-motion-of-a-cannonball-with-varying-drag

评论： Michael Madelaire 2018-12-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 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Michael Madelaire 2018-12-30

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/437689-projectile-motion-of-a-cannonball-with-varying-drag#answer_354431

在 MATLAB Online 中打开

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)')