物理学毕业论文设计运用MATLAB分析机械振动.doc

本科毕业论文题目: 运用Matlab分析机械振动 学院： 物理与电子科学学院 班级： 2008级物理二班 姓名： 指导教师： 职称： 实验师 完成日期： 2012 年 5 月 18 日运用Matlab分析机械振动作者：李钊 指导老师：陈翠萍（山西大同大学物理与电子科学学院，山西大同037009）摘要: 日常生活中所说的振动是一种周期性的运动。所谓周期性运动是指在时间上具有重复性或往复性的一种运动，是遍及自然界及社会科学界的一种运动方式。在物理学中，广义地说凡是描述物质运动的物理量，在某一数值附近做周期性的变化，都叫做振动.本文主要是例举了关于振动的典型实例，用Matlab语言编制计算机程序进行仿真以达到研究简谐振动以及振动的合成，振动的能量等目的，并在文章的最后简要地介绍了共振的危害与应用。关键字: Matlab语言 演示 振动 周期性 频率 合成 能量 共振 混沌现象 目 录一振动的概念与分类11.1狭义的振动11.1.1简谐振动的概念11.1.2简谐振动的动力学特征11.2广义的振动及振动的分类2二不同类型的振动合成及运用Matlab模拟演示32.1振动方向相同，频率相同的简谐振动的合成32.2 Matlab模拟振动方向相同，频率略有差异振动合成的“拍”现象42.3二个振动方向互相垂直的简谐振动的合成7三运用Matlab演示典型实例分析简谐振动的能量转换73.1简谐振动的系统机械能73.2 小球单摆的能量分析及Matlab演示：73.3Matlab演示弹簧振子：83.4简谐振动的能量曲线12四阻尼振动、受迫振动及位移共振154.1阻尼振动154.2受迫振动184.3位移共振21五、“不守规矩”的摆·混沌行为：225.1什么是“混沌”现象225.2依赖初值的两种情况22六、振动的危害226.1生产中接触到的振动源226.2振动引起的共振在历史事件中引起的危害226.3振动对人体器官的影响23七、共振创造了世界24致谢26一振动的概念与分类1.1狭义的振动1.1.1简谐振动的概念狭义的振动指的是机械振动，即力学系统中的振动。 振动(或机械振动)指的是物体在平衡位置附近往复运动。（琴弦、钟鼓、机械钟表的摆轮、发动机座、高耸的烟囱和固体晶格点阵的分子和原子都在振动） 振动是以波的形式传播的，机械振动的传播即机械波。 简谐振动是指质点在线性回复力作用下围绕平衡位置的运动。 1.1.2简谐振动的动力学特征（1）简谐振动过程中应该掌握的一些基本概念：振幅A：简谐运动物体离开平衡位置最大位移的绝对值A。振动的周期T：物体做一次完全振动所经历的时间。频率f：单位时间内物体所作的完全振动的次数。圆频率:一秒钟对应的圆心角。一次全振动对应的圆心角就是2（即360度）。相位：当振幅和频率一定时决定振动物体在任意时刻相对平衡位置的位移和速度的物理量。相位概念的重要性还在于比较两个简谐振动之间在“步调”上的差异。设有两个同频率的谐振动，它们的振动表达式为： 它们的相位差为： （2）以弹簧振子为例描述简谐振动的特点：弹簧自由伸展时，滑块的位置为原点O（即平衡位置）， 表示位移： 平衡位置O点线XXXXV图1-1-2 弹簧阵子由牛顿第二定律知： 由上式易知简谐振动的动力学方程 : 其解的形式为: 1.2广义的振动及振动的分类 从广义上说振动是指描述系统状态的参量（如位移、电压）在其基准值上下交替变化的过程。电磁振动习惯上称为振荡。力学系统能维持振动，必须具有弹性和惯性：由于弹性，系统偏离其平衡位置时，会产生回复力，促使系统返回来位置； 由于惯性，系统在返回平衡位置的过程中积累了动能，从而使系统越过平衡位置向另一侧运动。正是由于弹性和惯性的相互影响，才造成系统的振动。（1）按系统运动自由度分，有单自由度系统振动（如钟摆的振动）和多自由度系统振动。有限多自由度系统与离散系统相对应，其振动由常微分方程描述；无限多自由度系统与连续系统（如杆、梁、板、壳等）相对应，其振动由偏微分方程描述。（2）方程中不显含时间的系统称自治系统；显含时间的称非自治系统。（3）按系统受力情况分，有自由振动、衰减振动和受迫振动。（4）按弹性力和阻尼力性质分，有线性振动和非线性振动。（5）振动还可分为确定性振动和随机振动，后者无确定性规律。（如车辆行进中的颠簸）二不同类型的振动合成及运用Matlab模拟演示2.1振动方向相同，频率相同的简谐振动的合成合简谐振动的分振动方程为： 振动矢量合成如图示：图2-1振动矢量合成 用旋转矢量合成 合振幅矢量 合振动保持原振动方向不变。 合振动方程 由此易知：一个质点同时参与两个振动方向相同、频率相同的简谐振动，合振动仍为简谐振动。Matlab模拟编程如下：%两个同方向同频率的简谐振动的合成clear %清除变量a1=input('请输入第一个振动的振幅:); %第一个振动的振幅%a1=0.03; %参考值phi1=input('请输入第一个振动的初相的度数:');%第一个振动的初相%phi1=0; %参考值phi1=phi1*pi/180; %化为弧度a2=input('请输入第二个振动的振幅:'); %第二个振动的振幅%a2=0.04; %参考值phi2=input('请输入第二个振动的初相的度数:'); %phi1=0;90; phi2=phi2*pi/180; wt=linspace(0,4*pi); x1=a1*cos(wt+phi1); x2=a2*cos(wt+phi2); x=x1+x2; figure plot(wt,x1,'-.',wt,x2,'-',wt,x,'LineWidth',2)%画振动曲线set(gca,'XTick',(0:8)*pi/2) %设置横坐标刻度grid on %加网格fs=16; %字体大小title('同一直线上简谐振动的合成','FontSize',fs)%显示标题 2.2 Matlab模拟振动方向相同，频率略有差异振动合成的“拍”现象二个振动方向相同、频率略有差异的简谐振动，其合振动不为简谐振动，产生“拍”现象拍频为（为两分振动频率）“拍”现象合振动图像如下所示：x1tx2tx3t图2-2“拍”现象振动合成 Matlab演示“拍”现象：%拍的形成clear %清除变量d=10; %分母%d=15; %分母t=0:0.01:60; %时间向量w1=pi/2; %第一个角频率dw=pi/d; %角频率之差w2=w1+dw; %第二个角频率x1=cos(w1*t); %第一个位移x2=cos(w2*t); %第二个位移x=x1+x2; %合位移figure %创建图形窗口subplot(3,1,1) %选择子图plot(t,x1,t,x2,'-','LineWidth',2) %画位移曲线grid on %加网格leg1='itxrm_1/itArm=cositomegarm_1itt'%第一个图例字符串leg2='itxrm_2/itArm=cositomegarm_2itt'%第二个图例字符串legend(leg1,leg2) %图例tit='(itomegarm_1=pi/2,Deltaitomegarm=pi/',num2str(d),')'%标题一部分fs=16; %字体大小title('拍的形成' tit,'FontSize',fs) %加标题%xx1=cos(dw*t/2); %调幅线xx1=cos(w2-w1)*t/2); %调幅线(同上)%xx2=cos(w1+dw/2)*t); %无调幅的振动线xx2=cos(w2+w1)*t/2); %无调幅的振动线(同上)subplot(3,1,2) %选择子图plot(t,xx2,t,xx1,'r-','LineWidth',2) %画曲线grid on %加网格leg1='cos(itomegarm_2' %图例字符串的第一部分leg2='itomegarm_1)ittrm/2' %图例字符串的第二部分legend(leg1,'+',leg2,leg1,'-',leg2)%图例subplot(3,1,3) %选择子图plot(t,x1+x2,t,2*xx1,'r-',t,-2*xx1,'m-','LineWidth',2)%画曲线grid on %加网格xlabel('ittrm/s','FontSize',fs) %标记横坐标ylabel('itx/Arm=itxrm_1/itAit+itxrm_2/itA','FontSize',fs)%标记纵坐标2.3二个振动方向互相垂直的简谐振动的合成（1）若二振动频率相同，合振动轨迹一般为一椭圆. （2）若二振动频率成整数比，合振动轨迹为规则的稳定的闭合曲线，称利萨如图但若不成整数比，轨迹为不闭合的复杂曲线.三运用Matlab演示典型实例分析简谐振动的能量转换3.1简谐振动的系统机械能弹簧振子或扭摆振动系统中线性回复力为弹性力（或力矩），它们是保守力（或力矩），所以简谐振动系统的总机械能守恒。简谐振动的总机械能是简谐振动的动能与势能之和现以单摆、弹簧振子为例讨论振动系统的动能和势能的变化。3.2 小球单摆的能量分析及Matlab演示： 单摆的周期： 最高点与最低点的高度差：最高点时动能为0，最低点时势能为0所以振动的能量为： 抽象的问题具体化，运用Matlab演示单摆如下：%制作动画%挂摆横梁plot(-0.2;0.2,0;0,'color','y','linestyle','-',.'linewidth',10);%画初始位置的单摆g=0.98; %重力加速度,可以调节摆的摆速l=1;theta0=pi/4;x0=l*sin(theta0);y0=(-1)*l*cos(theta0);axis(-0.75,0.75,-1.25,0);axis('off'); %不显示坐标轴%创建摆锤head=line(x0,y0,'color','r','linestyle','.',.'erasemode','xor','markersize',40);%创建摆杆body=line(0;x0,0;y0,'color','b','linestyle','-',.'erasemode','xor');%摆的运动t=0;dt=0.01;while 1t=t+dt;theta=theta0*cos(sqrt(g/l)*t);x=l*sin(theta);y=(-1)*l*cos(theta);set(head,'xdata',x,'ydata',y);set(body,'xdata',0;x,'ydata',0;y);drawnow;end3.3Matlab演示弹簧振子：弹簧振子的弹性势能： 弹簧振子的动能： 弹簧振子的总机械能： 因为 ,均较易进行计算，所以计算动能时常用综上所述易知：任何简谐振动系统的机械能均可用下式计算简谐振动过程中，系统机械能守恒，但动能和弹性势能相互转化。简谐振动的振幅与机械能的关系。 理想弹簧振子的简谐振动Matlab编程如下:%理想弹簧阵子简谐运动%Clearrectangle('position',12,8.5,2,0.3,'FaceColor',0.5,0.3,0.4);axis(0,15,-1,10);%画顶板hold onplot(13,13,7,8.5,'r','linewidth',2);%画直线y=2:.2:7;M=length(y);x=12+mod(1:M,2)*2;x(1)=13;x(end-3:end)=13;D=plot(x,y); %弹簧C=0:.1:2*pi;r=0.35;t1=r*sin(C); F1=fill(13+r*cos(C),2+t1,'r');% 球set(gca,'ytick',0:2:9);set(gca,'yticklabels',num2str(-1:3');plot(0,15,3.3,3.3,'black');H1=plot(0,13,3.3,3.3,'y');% 句柄黄线Q=plot(0,3.8,'color','r');% 运动曲线;td=;yd=;T=0;text(2,9,'理想中的弹簧振子简谐振动','fontsize',16);set(gcf,'doublebuffer','on');while T<12;pause(0.2);Dy=(3/2-1/2*sin(pi*T)*1/2;Y=-(y-2)*Dy+7;Yf=Y(end)+t1;td=td,T;yd=yd,Y(end);set(D,'ydata',Y);set(F1,'ydata',Yf,'facecolor',rand(1,3);set(H1,'xdata',T,13,'ydata',Y(end),Y(end);set(Q,'xdata',td,'ydata',yd) ;T=T+0.1;End3.4简谐振动的能量曲线 能量曲线： 总机械能： 弹性势能能： 动能： 能量曲线如下图所示： Ex 图3-4简谐振动的能量曲线弹性势能与动能的平均值： 简谐振动中势能与动能的平均值相等且等于总机械能的一半。以弹簧振子为例运用Matlab模拟能量曲线：%弹簧振子的动能,势能和机械能曲线clear %清除变量n=4; %周期的个数t=linspace(0,2*pi)*n; %时间向量x=cos(t); %振子位置v=-sin(t); %速度ek=v.2; %动能ep=x.2; %势能figure %建立图形窗口subplot(2,1,1) %子图plot(t,x,t,v,'-','LineWidth',2) %画位移和速度曲线grid on %加网格axis tight %紧贴坐标轴fs=16; %字体大小title('简谐振动的位移和速度','FontSize',fs)%显示标题xlabel('itomegat','FontSize',fs) %标记横轴legend('位移itx/A','速度itv/omegaA')%图例set(gca,'XTick',(0:2*n)*pi) %设置横坐标刻度线text(0,0,'itomegarm=(itk/mrm)1/2','FontSize',fs)%显示角频率subplot(2,1,2) %子图plot(t,ek,'-',t,ep,'-.',t,ek+ep,'LineWidth',2)%画能量曲线grid on %加网格axis tight %紧贴坐标轴title('简谐振动的能量','FontSize',fs) %显示标题xlabel('itomegat','FontSize',fs) %标记横轴legend('动能itT/Erm_0','势能itV/Erm_0','机械能itE/Erm_0')%图例text(0,0.5,'itErm_0=itkArm2/2','FontSize',fs)%显示能量单位set(gca,'XTick',(0:2*n)*pi) %设置横坐标刻度线四阻尼振动、受迫振动及位移共振4.1阻尼振动以上讨论均假设质点或刚体的振动不受任何阻力，由于能量守恒，它们将永远振动下去。然后事实上，振动系统都受阻力作用，如无外界能量补偿，振动幅将不断减小而归于静止。振动系统因受阻力作振幅减小的运动，叫做阻尼振动。设 阻力 由牛顿第二定律得 令 由此可得其动力学方程 根据阻尼因数之不同，可将此方程解出三种可能的运动状态： x（1）欠阻尼状态 t 得质点的解 图4-1-1欠阻尼状态x（2）临界阻尼 t得质点的解 图4-1-2临界阻尼状态xt（3）过阻尼状态 得质点的解 图4-1-3过阻尼状态运用Matlab语言模拟阻尼运动：%阻尼运动的类型clear %清除变量t=0:0.25:20; %固有角频率与时间的乘积w0t向量(约化时间向量)%b=0:0.5:1.5; %阻尼因子与固有角频率的倍数向量(约化阻尼因子向量)b=0:0.25:1.25; %阻尼因子与固有角频率的倍数向量(约化阻尼因子向量)n=length(b); %曲线条数b(b=1)=1+eps; %将1值改为1加小量B,T=meshgrid(b,t); %约化阻尼因子和约化时间矩阵W=sqrt(1-B.2); %准角频率矩阵X=exp(-B.*T).*(cos(W.*T)+B./W.*sin(W.*T);%位移函数矩阵V=-exp(-B.*T).*sin(W.*T)./W; %速度函数矩阵%A=sqrt(B.2-1); %参数矩阵%X=exp(-B.*T).*(A+B).*exp(A.*T)+(A-B).*exp(-A.*T)/2./A;%位移函数矩阵(效果相同)%X=exp(-B.*T).*(cosh(A.*T)+B./A.*sinh(A.*T);%位移函数矩阵(效果相同)%V=-exp(-B.*T).*sinh(A.*T)./A; %速度函数矩阵(效果相同)f1=figure; %创建图形窗口%plot(t,X,'LineWidth',2) %画位移曲线族plot(t,X(:,1),'o-',t,X(:,2),'d-',t,X(:,3),'s-',t,X(:,4),'p-',. t,X(:,5),'h-',t,X(:,6),'<-') %画位移曲线族fs=16; %字体大小xlabel('itomegarm_0itt','FontSize',fs)%标记横坐标ylabel('itx/A','FontSize',fs) %标记纵坐标title('质点在不同阻尼下的运动曲线','FontSize',fs)%标题legend(repmat('itbeta/omegarm_0:',n,1),num2str(b')%加图例txt='itbeta/omegarm_0 小于1为欠阻尼,等于1为临界阻尼,大于1为过阻尼'%文本text(0,-0.7,txt,'FontSize',fs) %显示文本grid on %加网格f2=figure; %创建图形窗口%plot(t,V,'LineWidth',2) %画速度曲线族plot(t,V(:,1),'o-',t,V(:,2),'d-',t,V(:,3),'s-',t,V(:,4),'p-',. t,V(:,5),'h-',t,V(:,6),'<-') %画位移曲线族xlabel('itomegarm_0itt','FontSize',fs)%标记横坐标ylabel('itv/omegarm_0itA','FontSize',fs)%标记纵坐标title('质点在不同阻尼下的速度曲线','FontSize',fs)%标题grid on %加网格legend(repmat('itbeta/omegarm_0:',n,1),num2str(b')%加图例pause %暂停,可取图X1=; %位移矩阵清空V1=; %速度矩阵清空X2=; %位移矩阵清空V2=; %速度矩阵清空for i=1:n %按曲线循环 tm,XV=ode45('P5_7fun',t,1;0,b(i);%计算位移和速度 X1=X1,XV(:,1); %连接位移矩阵 V1=V1,XV(:,2); %连接速度矩阵 s='D2x+',num2str(2*b(i),'*Dx+x'%微分方程字符串sx=dsolve(s,'x(0)=1','Dx(0)=0'); %微分方程的符号解sv=diff(sx); %求速度的符号解x=subs(sx,'t',t); %位移v=subs(sv,'t',t); %速度 X2=X2,x' %连接位移矩阵 V2=V2,v' %连接速度矩阵end %结束循环figure(f1) %重开图形窗口hold on %保持图像plot(t,X1,'.',t,X2,'*') %画位移曲线figure(f2) %重开图形窗口hold on %保持图像plot(t,V1,'.',t,V2,'*') %画速度曲线%阻尼运动的二阶微分方程的函数function f=fun(t,x,flag,b)f= x(2); %速度 -2*b*x(2)-x(1); %加速度4.2受迫振动设质点受到三种力：弹性力-kx，阻尼力，周期性外力，亦称驱动力根据牛顿第二定律得，受迫振动的动力学方程： 令 ，得 解方程得： 开始时，受迫振动的振幅较小，经过一定时间后，阻尼振动即可忽略不计，质点进行由上式第二项所决定的与驱动力同频率的振动，称受迫振动的稳定振动状态，可表示如下：将此式带入方程得 ，xot图4-2 受迫振动自暂态发展为稳定振动.本图所示初始条件为t=0,初始条件影响暂态过程，不影响稳态振动Matlab编程演示受迫振动如下：%物体在平衡点从静止开始的受迫振动曲线(用解析式)clear %清除变量b=input('请输入约化阻尼因子(01):'); %键盘输入约化阻尼因子%b=0.1; %参考值if b<=0|b>=1 return,end %不符合条件则不向下执行程序w=sqrt(1-b2); %约化阻尼圆频率s='请输入约化驱动力圆频率(约化阻尼圆频率为',num2str(w),'):'%提示字符串W=input(s); %键盘输入约化驱动力圆频率%W=2;6;1;0.6; %参考值if W=1 W=1-eps;end %如果为1则改小一点tm=30; %最大时间t=0:0.001:tm; %时间向量a1=sqrt(w2*(W2-1)2+b2*(W2+1)2)/w/(W2-1)2+4*b2*W2);%阻尼振幅phi=atan2(b*(W2+1),w*(W2-1); %阻尼振动初相a2=1/sqrt(W2-1)2+4*b2*W2); %等幅振动振幅PHI=atan2(-2*b*W,1-W2); %等幅振动初相x1=a1*exp(-b*t).*cos(w*t+phi); %阻尼振动函数x2=a2*cos(W*t+PHI); %等幅振动函数x=x1+x2; %合成振动xm=max(abs(x); %最大值%-figure %创建图形窗口subplot(3,1,1) %选子图plot(t,x1,'LineWidth',2) %画曲线grid on %加网格axis(0,tm,-xm,xm) %设置曲线范围fs=12; %字体大小title('减幅振动的位移时间曲线','FontSize',fs)%标题ylabel('itxrm_1/itArm_0','FontSize',fs)%标记纵坐标txt='itbeta/omegarm_0=',num2str(b);%阻尼因子字符串txt=txt,',itomega/omegarm_0=',num2str(w);%连接阻尼圆频率text(0,xm,txt,'FontSize',fs) %标记阻尼因子和阻尼圆频率subplot(3,1,2) %选子图plot(t,x2,'LineWidth',2) %画曲线grid on %加网格axis(0,tm,-xm,xm) %设置曲线范围title('等幅振动的位移时间曲线','FontSize',fs)%标题ylabel('itxrm_2/itArm_0','FontSize',fs)%标记纵坐标txt='itOmega/omegarm_0=',num2str(W);%驱动力圆频率字符串text(0,xm,txt,'FontSize',fs) %标记驱动力圆频率subplot(3,1,3) %选子图plot(t,x,'LineWidth',2) %画曲线grid on %加网格axis(