MATLAB 在存在空气阻力情况下自由落体的模型

 模拟的是水滴模型，所以在空气阻力方面我是用的1/2c*p*s*v，是使用的v的一次方，如果是高速的话应该用二次或者三次方，利用不定积分与微分方程，将速度与t的关系式表达出来，然后我使用了微分法，dt应该可以更小一点会更准确，不断将v*dt相加，然后画图

clc;                                                        
clear all;
p=1;
r1=0.001;
r2=input('请输入第二个半径是第一个的几倍，r2=\n')/1000;
V1=(4*pi*r1^3)/3;
V2=(4*pi*r2^3)/3;
s1=pi*r1^2;
s2=pi*r2^2;
k=0.00018102;
T=5;
g=9.8;
H0=1/2*g*T^2;
H1=0;                                                      
H2=0;
dt=0.1;
syms v
for t=0:dt:T
    answ1=solve(p*r1/k*log(g/(g-k*v/(p*r1))) == t,v);
    H1=H1-answ1*dt;
    answ2=solve(p*r2/k*log(g/(g-k*v/(p*r2))) == t,v);
    H2=H2-answ2*dt;
    subplot(1,2,1);
    plot(0,H1,'.r','markersize',20);
    ylabel('height(m)');
    ylim([-H0,0]);                                          
    subplot(1,2,2);
    plot(0,H2,'.r','markersize',20000*r2);
    ylabel('height(m)');
    ylim([-H0,0]);
    drawnow;
end

空气阻力系数可能不太准确，需要修改一下。