《MATLAB计算单摆的周期和动画.docx》由会员分享，可在线阅读，更多相关《MATLAB计算单摆的周期和动画.docx（13页珍藏版）》请在三一办公上搜索。
1、MATLAB计算单摆的周期和动画单摆振动的周期和运动规律 问题 单摆振动的周期和运动规律 (1)求单摆的周期与角振幅的关系。 (2)演示单摆的振动。 数学模型 (1)如A5.1图所示，设摆锤质量为m，角位置为，摆锤的运动方程为 d2qml2=-mgsinq， dtd2qg=-sinq， (5.1.1) 即 2dtl在小角度的情况下，sin ，可得 O l T ft mg d2q2+w0q=0， (5.1.2) 2dt其中w0=g/l，0为圆频率。可知：单摆在小角度时作简谐振动，A5.1图 小角度周期为 T0=2w0=2l。 (5.1.3) g可见：在小角振动的情况下，单摆的周期与角振幅无关，这
2、称为单摆的等时性。 摆锤的角速度为 = d/dt，因此 d2qdwdwdqdw=w， 2dtdtdqdtdq由(5.1.1)式可得 gwdw=-sinqdq， l积分得 12gw=cosq+C， 2l当t = 0时， = 0， = m，可得C = -gcosm/l。因此角速度大小为 dq2g=w=(cosq-cosqm)。 (5.1.4) dtl注意：角速度是单位时间内角度的变化率d/dt，圆频率是简谐运动中2时间内周期性运动的次数2/T，它们常用字母表示，单位也相同，但意义不同。单摆的周期为 1 lT=42gqm0dq2=T0cosq-cosqmqm0dq。 (5.1.5) cosq-cos
3、qm对于任何角振幅m，通过数值积分和符号积都能计算周期。 利用半角公式可得 1T=T0设 qm0dqsinqm/2-sinq/222。 k=sinqm， (5.1.6) 2并设ksinx = sin(/2)，因此kcosxdx=1qcosdq，可得 22/21T=T0/202=T0cos(q/2)k2-k2sin2x/22kcosxdx0dx1-sinq/22， 2即 T=T00/2dx1-ksinxdx1-ksinx2222。 (5.1.7) 这是椭圆积分。第一类完全椭圆积分定义为 K(k)=周期为 0， (5.1.8) T=T02K(k)。 (5.1.9) 算法对于任何一个角振幅m，利用(
4、5.1.5)式，通过MATLAB数值积分指令quadl和符号积分指令int都可计算单摆的周期。利用MATLAB的完全椭圆积分指令ellipke也能计算单摆的周期。不过，MATLAB定义的一类完全椭圆积分定义为 /2K(m)=周期可表示为 0dx1-msinx2， (5.1.8*) T=T0其中 2K(m)。 (5.1.9*) m=sin2程序zqy5_1circle.m如下。 2 qm。 (5.1.6*) 2%单摆的周期(用内线函数数值积分和符号积分与椭圆积分比较) clear %清除变量 theta=1:179; %角振幅向量 thm=theta*pi/180; %角度化为弧度数 T0=el
5、lipke(sin(thm/2).2)*2/pi; %用椭圆积分计算单摆精确周期 T1=; %周期向量置空 T2=; %周期向量置空 n=length(thm); %角度个数 j=1:5:n; %下标向量(每隔5度取一个角度) for i=j %按下标循环 s=1./sqrt(cos(x)-cos( num2str(thm(i) );%被积函数字符串 f=inline(s); %被积的内线函数 t=quadl(f,0,thm(i)*sqrt(2)/pi; %用数值积分求周期 T1=T1,t; %连接周期 t=double(int(s,0,thm(i)*sqrt(2)/pi;%用符号积分求周期
6、T2=T2,t; %连接周期 end %结束循环 figure %创建图形窗口 plot(theta,T0,theta(j),T1,.,theta(j),T2,o)%画周期曲线 grid on %加网格 fs=16; %字体大小 xlabel(itthetarm_m/circ,FontSize,fs)%x标签 ylabel(itT/Trm_0,fontsize,fs) %y标签 title(单摆的周期与角振幅的关系,FontSize,fs)%标题 legend(椭圆积分,数值积分,符号积分,2)%图例 text(0,2,itTrm_0=2pi(itl/grm)1/2,FontSize,fs)%
7、小角单摆周期 图示如P5_1图所示，数值积分和符号积分与椭圆积分的结果完全吻合，而椭圆积分的效率最高，应用也最简单。当角振幅在20以内时，单摆的周期几乎不变。当角振幅在20到40之间时，单摆的周期稍有增加。当角振幅在大于40时，单摆的周期显著增加。当角振幅接近180时，单摆的周期急剧增加。 解析(2)为了演示单摆的振动，需要求微分方程中角度的数值解。摆锤的坐标为 x = lsin，y = -lcos， (5.1.9) 其中，x轴取向右的方向为正，y轴取向上的方向为正。 P5_1图 单摆的振动可与简谐振动进行比较。简谐振动的角位移可用余弦函数表示 = mcost。 (5.1.10) 简谐振动的角
8、速度为 w=dq=-qmwsinwt。 (5.1.11) dt 3 算法取t* = 0t为量纲时间，(5.1.1)式可化为 d2q=-sinq， (5.1.1*) dt*2取x(1) = ，x(2) = d/dt*，可得 dx(1)dx(2)=x(2)=-sinx(1)。 (5.1.1*) ，*dtdt初始条件是x(1) = m，x(2) = 0。由于 x(2)=dx(1)dt*=dqwdt=ww， 0因此x(2)表示以0为单位的角速度。由于 t=w0tw=w0t2T0， 0因此，无量纲时间除以2就是以小角周期T0为单位的时间。 取l为坐标单位，则摆锤的坐标可表示为 x* = x/l = si
9、n，y* = y/l = -cos。 简谐振动的角位移可表示为 q=q2tmcos(1T/T)=q1mcos(T*w0t)。 0T0其中T* = T/T0是约化周期。简谐振动的角速度可表示为 w=-q21Tsin(12tT)=-wqm1m0*sin(*w0t)。 0T/T0T/T00TT 程序zyq5_2ode.m如下。 %单摆振动的动画 clear %清除变量 thetam=input(请输入单摆角振幅的度数:);%键盘输入角振幅 thm=thetam*pi/180; %化为弧度 T=ellipke(sin(thm/2).2)*2/pi; %用椭圆积分计算精确周期 wt=linspace(0
10、,2*pi*T); %无量纲时间向量 wt,X=ode45(zqy5_2fun,wt,thm,0);%计算角度和角速度 th=X(:,1); %角度向量 w=X(:,2); %角速度向量 figure %创建图形窗口 subplot(2,1,1) %选子图 plot(wt/2/pi,th*180/pi,wt/2/pi,thetam*cos(wt/T),.)%画角度曲线 grid on %加网格 axis tight %曲线贴框 fs=16; %字体大小 title(单摆的角度,FontSize,fs) %标题 4 (5.1.9*) (5.1.10*) (5.1.11*) xlabel(itt/
11、Trm_0,FontSize,fs) %横坐标 ylabel(角度itthetarm/circ,FontSize,fs)%纵坐标 text(0,0,itTrm_0=,num2str(T),FontSize,fs)%显示周期 subplot(2,1,2) %选子图 plot(wt/2/pi,w,wt/2/pi,-thm/T*sin(wt/T),.)%画角速度曲线 grid on %加网格 axis tight %曲线贴框 title(单摆的角速度,FontSize,fs) %标题 xlabel(itt/Trm_0,FontSize,fs) %横坐标 ylabel(角速度itomega/omega
12、rm_0,FontSize,fs)%纵坐标 pause figure %建立图形窗口并取句柄 axis equal off %不显示坐标 title(单摆的振动,FontSize,20) %显示标题 hold on %保持图像 plot(0,0,o) %画悬点 plot(0;0,0;-1.05,-.) %画竖虚线 plot(exp(i*(linspace(-thm,thm)-pi/2),m-)%画轨迹 x=sin(thm); %摆球的起始横坐标 y=-cos(thm); %摆球的起始纵坐标 plot(0 -x,0 y,r:) %画连线 %pole=plot(0 x,0 y,r,LineWidt
13、h,2,EraseMode,xor);%取摆线的句柄 pole=plot(0 x,0 y,r,LineWidth,2);%取摆线的句柄 %ball=plot(x,y,c.,MarkerSize,50,EraseMode,xor);%取摆球的句柄 ball=plot(x,y,c.,MarkerSize,50); %取摆球的句柄 pause %暂停 while 1 %无限循环 for i=1:length(wt) %按时间循环 x=sin(th(i); %摆球的横坐标 y=-cos(th(i); %摆球的纵坐标 set(ball,XData,x,YData,y) %设置球的坐标 set(pole,
14、XData,0 x,YData,0 y)%设置摆绳的坐标 drawnow %刷新屏幕 pause(0.001) %延时 end %结束循环 if get(gcf,CurrentCharacter)=char(27) break;end%按ESC键则中断循环 end %结束循环 程序要调用函数zqy5_2fun.m。 %单摆运动规律的函数文件 function f=fun(t,x) f= x(2); %角速度表达式 -sin(x(1); %角加速度表达式 说明程序执行时从键盘输入角振幅，显示单摆之后按任一键，就演示单摆的振动，直到按ESC键为止。 5 图示(1)如果角振幅为60，如P5_2a图所示，尽管角振幅比较大，但是单摆的周期与小角单摆的周期相比，增加得很少。单摆的角位移与简谐振动相近。单摆的初始状态如P5_2b图所示。 (2)如果角振幅为120，单摆振动的角位移和角速度如P5_2c图所示。这时大角单摆的周期显著大于小角单摆的周期。单摆的初始状态如P5_2d图所示。 (3)如果角振幅只有5或小于5，就成了小角单摆(图略)。 P5_2a图 P5_2c图 P5_2b图 P5_2d图 6
