MATLAB计算平抛阻尼运动.doc

小球在空气中平抛运动的轨迹问题一小球在空气中作平抛运动，初速度为v0，所受的阻力与速率成正比：f = -kv，k称为阻力系数。画出小球运动的轨迹。 数学模型mgxxOv0fB7.1图yx如B7.1图所示，小球受到重力mg，方向向下；空气阻力f，方向与速度方向相反。根据牛顿第二定律可列出直角坐标方程，。 (7.1.1)由于vx = dx/dt，vy = dy/dt，上式可化为，。 (7.1.2)分离变量得，积分得，当t = 0时，vx = v0，vy = 0，可得Cx = -v0，Cy = -g，因此，。 (7.1.3)当t = 0时，x = 0，y = 0，积分上式可得，。 (7.1.4)这是小球的运动方程，也是以时间t为参数的轨道方程。算法方法一：用解析式。取t0 = m/k为时间单位，取s0 = m2g/k2为坐标单位，坐标方程可表示为x = s0v0*1 exp(-t*)，y = s0t* + exp(-t*) - 1。 (7.1.4*)其中，t* = t/t0，为约化时间；v* = kv0/mg，是小球的约化初速度。程序zqy6_3ode.m的第一部分如下。%阻力与速度成正比的平抛运动的轨迹(用解析解)clear %清除变量v0=input(请输入水平初速度kv0/mg:); %键盘输入水平初速度t=0:0.1:4; %时间向量x=v0*(1-exp(-t); %有空气阻力的横坐标y=t+exp(-t)-1; %有空气阻力的纵坐标fs=16; %字体大小figure %开创图形窗口plot(x,y) %画轨迹grid on %加网格xlabel(itx/srm_0,FontSize,fs) %标记横坐标ylabel(ity/srm_0,FontSize,fs) %标记纵坐标title(阻力与速度成正比的平抛运动的轨迹,FontSize,fs)%标题text(0,max(y),itkvrm_0/itmgrm=,num2str(v0),FontSize,fs)%约化初速度axis ij equal %原点设在左上角并使坐标间隔相等说明程序执行时，从键盘输入约化初速度，例如1等值。方法二：用两个微分方程的数值解。微分方程组(7.1.1)式可化为，。其中，t0 = m/k，V0 = gt0 = mg/k。取t* = t/t0，vx* = vx/V0，vy* = vy/V0，可得，。 (7.1.1*)取v(1) = vx*，v(2) = vy*，可得两个一阶方程组，。在初始时刻，小球约化初速度为，而v(2) = 0。在任意时刻，v(1)和v(2)表示约化速度。根据速度值，利用指令cumtrapz指令可求坐标。程序zqy6_3ode.m的第二部分如下。%阻力与速度成正比的平抛运动的轨迹(求两个一阶微分方程的数值解)t0,V=ode45(zqy6_3fun1,t,v0,0);%求微分方程的数值解x=cumtrapz(V(:,1)*t(2); %横坐标y=cumtrapz(V(:,2)*t(2); %纵坐标hold on %保持图像plot(x,y,r.) %画轨迹程序在执行时将调用一个函数zqy6_3fun1.m。%阻力与速度成正比的平抛运动的加速度函数function f=fun(t,v) f=-v(1); %水平加速度 1-v(2); %竖直加速度说明用两个常微分方程的数值解可求出速度，还需要求位置坐标才能画轨迹。方法三：用四个微分方程的数值解。微分方程组(7.1.1)式可化为，。取t* = t/t0，x* = x/s0，y* = y/s0，由于t0 = m/k，s0 = m2g/k2，可得，。 (7.1.1*)取r(1) = x*，r(2) = y*，r(3) = dx*/dt*，r(4) = dy*/dt*，可得四个一阶方程组； ，；。在初始时刻，小球约化位移为r(1) = r(2) = 0，初始约化速度为，r(4) = 0。在任意时刻，r(1)和r(2)表示约化坐标，r(3)和r(4)表示约化速度。程序zqy6_3ode.m的第三部分如下。%阻力与速度成正比的平抛运动的轨迹(求四个一阶微分方程的数值解)t0,R=ode45(zqy6_3fun2,t,0,0,v0,0);%求微分方程的数值解plot(R(:,1),R(:,2),ko) %画轨迹程序在执行时将调用一个函数zqy6_3fun2.m。%阻力与速度成正比的平抛运动的函数function f=fun(t,r) f= r(3); %水平速度 r(4); %竖直速度 -r(3); %水平加速度 1-r(4); %竖直加速度说明物体在平面上运动时，运动的微分方程组一般可化四个一阶常微分方程，要注意速度和加速度的排列顺序。注意如果取r(1) = x*，r(2) = dx*/dt*，r(3) = y*，r(4) = dy*/dt*，四个一阶方程组需要修改如下； ，；。在初始时刻，小球约化位移为r(1) = r(3) = 0，初始约化速度为，r(4) = 0。在任意时刻，r(1)和r(3)表示约化坐标，r(2)和r(4)表示约化速度。函数文件需要修改如下(zqy6_3fun3.m)%阻力与速度成正比的平抛运动的函数function f=fun(t,r) f= r(2); %水平速度 -r(2); %水平加速度 r(4); %竖直速度 1-r(4); %竖直加速度主程序需要修改如下t0,R=ode45(zqy6_3fun3,t,0,v0,0,0);%求微分方程的数值解plot(R(:,1),R(:,3),ko) %画轨迹方法四：用微分方程的符号解。公式(7.1.1*)可化为二阶微分方程，。 (7.1.1*)据此可求解微分方程的符号解。程序zqy6_3ode.m的第四部分如下。%阻力与速度成正比的平抛运动的轨迹(用微分方程的符号解)sx,sy=dsolve(D2x+Dx,D2y-1+Dy,Dx(0)=v0,Dy(0)=0,x(0)=0,y(0)=0)%求微分方程的符号解x=subs(sx,v0,v0); %替换初速度x=subs(x,t,t); %替换时间形成横坐标y=subs(sy,v0,v0); %替换初速度y=subs(y,t,t); %替换时间形成纵坐标plot(x,y,ms) %画轨迹说明符号解有两个结果sx =v0-v0*exp(-t)sy =exp(-t)+t-1这