1、第七章 动力学与振动,7.1 轨迹 7.2 单自由度系统 7.3 多自由度系统,第7章 振动,7.1 轨迹,举例说明:重力场中有两个物体,其中质量为m2的物体固定,而质量为m1的物体绕m2做平面圆周运动。做圆周运动的m1物体的轨道半径用变量r表示,角度用变量表示。,两物体系统,卫星绕地球转动时，m2等于地球的质量，m1等于卫星的质量，r为卫星球心与地球球心间的距离。其运动轨迹由下列方程组决定：,引入状态变量:,第7章 振动,建立函数文件orbit.m function xd=orbit(t,x) xd=x(2);x(1)*x(4)2-4.0*pi2/x(1)2;x(4);-2.0*x(2)*x
2、(4)/x(1); 三组初始条件(t=0)：,由初始条件建立执行文件menu71.m initcond=2 0 0 1.5;1 0 0 2*pi;2 0 0 4; tspan=linspace(0,5,1000); options=odeset(RelTol,1e-6,AbsTol,1e-6 1e-6 1e-6 1e-6); lintype=-. -. -.; for i=1:3t,x=ode45(orbit,tspan,initcond(i,:),options);polar(x(:,3),x(:,1),lintype(2*(i-1)+1:2*i);hold on end text(0.5,
3、-1.2,椭圆轨迹); text(-1.2,1,圆轨迹); text(1.75,2,双曲线轨迹);,第7章 振动,运行结果,第7章 振动,第7章 振动,7.2 单自由度系统,一、概述,1、力学模型,其中：振体质量为m，弹簧的线性系数为k，非线性系数为， 阻尼系数为c，外力F(t）。,弹簧质量阻尼系统,2、运动微分方程,第7章 振动,用x表示系统的位移，则运动微分方程为：,引入新变量转化状态空间方程形式：,二、线性系统的自由振动,1、运动微分方程,当线性项为0时，得到线性振动系统的自由振动方程。,2、MATLAB求解,对应的函数文件FreeOcillation.m function xdot=F
4、reeOcillation(t,x,dummy,zeta) xdot=x(2);-2.0*zeta*x(2)-x(1);,第7章 振动,三种阻尼系数（1）阻尼系数为0.1时是欠阻尼情况；（2）阻尼系数为1时是临界阻尼情况；（3）阻尼系数为5时是过阻尼情况。,由初始条件（位移和速度均为1时）建立执行文件menu72.m zeta=0.1 1.0 5.0; tspan=linspace(0,40,400); lintype=-b -r -r; for i=1:3t,x=ode45(FreeOcillation,tspan,1 1,zeta(i);subplot(2,1,1); plot(t,x(:
5、,1),lintype(2*(i-1)+1:2*i);hold onsubplot(2,1,2); plot(x(:,1),x(:,2),lintype(2*(i-1)+1:2*i);hold on end subplot(2,1,1);,第7章 振动,xlabel(Time( tau); ylabel(Displacement x( tau); title(Displacement as a function of( tau); axis(0 40 -2.0 2.0); text(2.7,-1.3,阻尼系数=0.1);text(3.6,-0.1,1.0); text(3.6,1.0,5.0)
6、; subplot(2,1,2); xlabel(Displacement); ylabel(Velocity); title(Phase portrait); axis(-2.0 2.0 -2.0 2.0); text(0.7,-1.25,阻尼系数=0.1);text(0.8,-0.65,1.0); text(0.8,0.1,5.0);,运行结果,第7章 振动,第7章 振动,三、线性系统的强迫振动,1、运动微分方程,2、MATLAB求解,对应的函数文件ForceOcillation.m function xdot=ForceOcillation(t,x,dummy,zeta,Omega,x0
7、) xdot=x(2);-2.0*zeta*x(2)-x(1)+x0*cos(Omega*t);,为了获得频谱图,建立函数文件AmplitudeSpectrum.m functionf,amplitude=AmplitudeSpectrum(yy,Fs,Nstart,N); f=(Fs*(0:N-1)/N)*2.0*pi; amplitude=abs(fft(yy(Nstart:Nstart+N),N)/N; 采样速率30/6000=0.005,则采样频率1/0.005=200,这个频率远远超出了必须达到的采样频率,结果显示截短频谱图,需设置Nstart=3200,N=211=2048。,t
8、= 0:0.001:0.6; x = sin(2*pi*50*t)+sin(2*pi*120*t); y = x + 2*randn(size(t); plot(y(1:50) title(Signal Corrupted with Zero-Mean Random Noise) xlabel(time (seconds),第7章 振动,fft的应用,Fft变换程序,Y = fft(y,512); Pyy = Y.* conj(Y) / 512; f = 1000*(0:256)/512; plot(f,Pyy(1:257) title(Frequency content of y) xlab
9、el(frequency (Hz),第7章 振动,时域运行结果,频域运行结果,编制执行文件menu72f.m zeta=0.4;Omega=3.0;x0=50; tspan=linspace(0,30,6000); options=odeset(RelTol,1e-8,AbsTol,1e-8); lintype=-b;t,x=ode45(ForceOcillation,tspan,0 0,options,zeta,Omega,x0);subplot(2,1,1);plot(t,x(:,1);axis(0 30 -8 8);hold onsubplot(2,1,2);,第7章 振动,yy=x(:
10、,1);N=2048;Nstart=3200;Fs=200; f,Amplitude=AmplitudeSpectrum(yy,Fs,Nstart,N);semilogy(f(1:40),2*Amplitude(1:40); xlabel(Frequency); ylabel(Amplitude); title(Response spectrum of a linear system);hold on subplot(2,1,1); xlabel(Time( tau); ylabel(Displacement x( tau); title(Response of a linear system
11、); hold on,第7章 振动,运行结果,第7章 振动,频率响应 阶跃响应 脉冲响应,Bode函数: g = tf(1 0.1 7.5,1 0.12 9 0 0); bode(g),第7章 振动,第7章 振动,bode(g,0.1 , 100),nyquist(g),impulse(g) 脉冲响应,step(g) 阶跃响应,一、建立系统的运动微分方程,1.用刚度影响系数法建立振动微分方程,可简写为,第7章 振动,7.3 自由振动模态及固有频率,例：三质量m1，m2，m3串联于弹簧k1，k2，k3上，试列出自由振动微分方程。,解：求刚度系数。给x1以单位位移，x2与x3保持不动，即x1=1，
12、x2=x3=0；要产生这样的位移状态，在各点加的力就是k11，k21，k31，显然有： k11=k1+k2，k21=k2，k31=0（得到第一列） 再令x2=1，x1=x3=0， 则有 k12=k2，k22=k2+k3，k32=k3,第7章 振动,于是得到刚度矩阵：,由于系统有对应于广义坐标x1，x2，x3的集中质量，故质量矩阵是对角阵，于是得到自由振动微分方程：,第7章 振动,2.用柔度影响系数法建立振动微分方程,柔度影响系数法就是力法，它所建立方程是各点的位移协调方程。柔度影响系数cij定义：在j点作用有单位力，而其他各点没有力，所引起i点的位移。,从右图看出，当j点作用单位力时，在梁上各
13、点产生的位移c1j，c2j，cij，cnj。我们暂时不管其他点，只研究i点引起位移cijPj。,同样，在1点加力P1，在i点产生位移ci1P1；,这就是i点的位移协调方程，用同样方法可以得到梁上各点的位移协调方程：,设在力系P1，P2Pn作用下i点的位移是yi，那么必有：,同样，在2点加力P2，在i点产生位移ci2P2；,同样，在j点加力Pj，在i点产生位移cijPj；,同样，在n点加力Pn，在i点产生位移cinPn；,第7章 振动,写成矩阵形式：,（1）,其中,称为柔度矩阵。由功的互等定理可知cij=cji，因此c是对称矩阵。我们前边讲过，自由振动时 ，即,第7章 振动,于是（1）式可写成：
14、,此即用柔度影响系数法建立的多自由度系统自由振动微分方程。,或简写成：,（2）,第7章 振动,讨论 ：,就是说，刚度矩阵和柔度矩阵是互为逆矩阵。,1、（2）式可写成,等式两边同乘 ，则,即,与前边得到的 比较，可见：,按逆阵性质：,第7章 振动,2、方程（2）与建立的方程形式不同，实质是一样的。,解：首先计算柔度系数，用莫尔法：,例：一悬臂梁，固定三个集中质量，梁的抗弯刚度EI为常数，梁质量不计，求：试用柔度系数法列微分方程。,(单位：长度/力),注：1、Mi在i加单位力的弯矩方程，Mj在j加单位力的弯矩方程；2、原公式中 号指弯矩方程不连续时，分段积分后求和。,第7章 振动,于是：,第7章
15、振动,故柔度矩阵为,即,自由振动微分方程：,第7章 振动,关于固有频率及固有振型的概念，我们已经学习过。研究多自由度系统的这些固有特征时, 采用矩阵为数学工具。,二、固有频率与固有振型,注意:对于一个正定系统（即系统没有刚体位移），用上述各种方法建立微分方程都可以，有时用柔度系统法更为方便；但对于半正定系统（系统有刚体位移），则柔度系统没有意义（加单位力后系统不能维持平衡而产生刚体运动），只能用形如My+ky=0那样的方程。,第7章 振动,对于有n个自由度q1，q2qn的系统，其自由振动微分方程的一般形式是：,或写成矩阵形式,这是一组二阶常系数、线性常微分方程。设系统偏离平衡位置作自由振动时，
16、存在各qi按同一频率，同一相位作简谐振动的特解：,（1）,第7章 振动,这是一组关于振幅A1，A2An的线性齐次代数方程组，要使A1，A2An有非零解（零解对应不振动），必有系数行列式为零，即：,代入方程组（1），消去公因子 后，得,（2）,或,第7章 振动,或,此即多自由度系统的特征方程（或称频率方程）,在这个方程中，kij，mij都是已知的，未知的只有，展开后是一个关于2的n次代数方程，从中可以求得2的n个根。这些根就是系统的固有频率。由小到大排列：1，2n叫做基阶，二阶n阶固有频率。,第7章 振动,把每一个固有频率例如代回方程组（2），就可以求得该阶固有振型：,这里解释一下。在二自由度系
17、统时，求得了固有频率1，2之后，把1（或2）代回振型方程组求该阶振型时，是代回两个方程中的任意一个，求得，及 。就是说，振型方程组的两个方程不独立，而是线性相关。对于多自由度系统，也是一样，因为振型方程组（2）的系统行列式为零，即 ，就说明矩阵 的秩rn，至少有一个方程与其他方程线性相关。我们在求固有振型时，应划出一个方程，例如最后一个，而把剩下的n-1个方程中某一项移到等式右边，例如把A1项移过去。那么解出的A2，A3An都包含A1，也就是说，求出了A2，A3An各自对A1的比值，那么它们之间的比值自然也就确定了。,第7章 振动,即该阶固有振型：,显然：,代回我们设的特解,注：上标（1）表示
18、第一阶振型,或 第1阶主振动,第7章 振动,第1阶振型：,同样，当把2，3n分别代入振型方程（2）时，可求得其他各阶固有振型。即：,就是说，当系统按第一阶固有频率作振动时，各点振幅 之间具有确定的比值，而且系统各点之间在任一瞬时也具有同样确定的比值，这说明系统有一定的振动形态，称为第一阶固有振型。,第2阶振型：,第n阶振型：,于是也可写出各阶主振动。,第7章 振动,方程（1）的通解是这n个主振动的叠加：,写成矩阵形式：,（3）,式中共有2n个特定常数 ：及 ，注：对于每一阶振型来说，各阶振幅之间比值已定，只要知道其中一个就可求得其他振幅，所以，每一振型中只有一个待定振幅。这些常数由初始条件（t
19、=0时的 及 ）来确定。,第7章 振动,由（3）式可见，具有n个自由度的无阻尼系统自由振动一般由n个不同频率的主振动组成。在每一主振动中，系统各点的位移始终保持一定的比例关系（等于各点振幅之比），即系统有一定的振动形态，称为主振型（或固有振型）。可以任取出其中一点的振幅（例如A1）等于1，使其标准化，这样的振型，称为标准化了的固有振型。,解：取三个扭转角1，2，3为广义坐标。,例一、悬臂轴抗扭刚度GJ，质量不计，其上固定三个圆盘，转动惯量是I。求：系统的固有频率，固有振型。,第7章 振动,代入拉氏方程,得系统自由振动微分方程：,第7章 振动,或写成矩阵形式：,把得到的M、k矩阵代入振型方程,注
20、：或设 ,设 ，则上式成为：,注：可以直接代入频率方程 ,第7章 振动,频率方程。,由一元三次方程求根公式解得：,展开为：,1=0.198，2=1.555，3=3.247,代入 式中，求得三个固有频率,第7章 振动,代回振型方程中的前二个（用），把A3移到方程右端：,将上面求得的1，2，3分别代入，就得到三个固有振型：,取A1=1（即对标准化），得,即：,第7章 振动,第7章 振动,编制m文件. k=2 -1 0;-1 2 -1;0 -1 1; m=1 0 0;0 1 0;0 0 1; vibrationmode,eigenvalues=eig(k,m) v1=vibrationmode(:,
21、1) v2=vibrationmode(:,2) v3=vibrationmode(:,3) b1=v1./v3 b2=v2./v3 b3=v3./v3 for i=1:3 subplot(3,1,i) plot(0 1 2 3,0 b3(i) b2(i) b1(i) hold on end,运行结果,第7章 振动,应用MATLAB求解,建立系统的运动微分方程,在前一节看到，为了求多自由度系统的固有频率和固有振型，必须解高次代数方程（频率方程或称特征方程），其次数与系统自由度数目一致。当自由度数增加时，解次数很高的代数方程，是很困难的。所以一般用近似的数值方法。这里介绍的是一种数值方法“矩阵迭
22、代法”。,多自由度相同自由振动微分方程一般形式：,设特解为：,第7章 振动,三、用矩阵迭代法求固有频率和固有振型,代入方程得振型方程,用k-1前乘各项，移项后得到,引入符号,动力矩阵，,就得到迭代格式（振型方程）：,对于一个给定的系统，动力矩阵D是可知的。矩阵迭代法就是用逐次迭代、逐次逼近的方法求解振型方程，以确定特征值 和特征向量（固有振型）A。,第7章 振动,集中质量法概念,第7章 振动,例:集中质量法求简支梁的模态,力学模型:有一简支梁，其上有一个附加质量（位置和尺寸均在图上标注），梁长L=2m ,截面尺寸b*h为0.04m*0.02m ,截面惯性矩为2.67* ，密度为7920kg/，
23、质量块的质量M为3.0kg 弹性模量E=210GPa,离散,第7章 振动,一.力学模型及集中质量法简化,质量矩阵,柔度矩阵,矩阵 = 称作系统的动力矩阵,=,第7章 振动,二、MATLAB计算,m=1 0 0;0 1.94 0;0 0 1; h=0.00000195*3.186; b=9 11 7;11 16 11;7 11 9; f=h.*b; a=f*m; d=1 1 1; for i=1:16; y=a*d; d=(1/y(3,1)*y; end v=d s=m*y;,第7章 振动,w1=1/sqrt(s(3,1) subplot(3,1,1) plot(0 2,0 0) hold on
24、 plot(0 0.5 1 1.5 2,0 v(1,1) v(2,1) v(3,1) 0) p=v*m*v; q=v*v*m; j=a-(s(3,1)/p)*q; d=1 1 -1; % x(1)=d; for i=1:16; y=j*d; % l=y(i) d=(1/y(3,1)*y; end,第7章 振动,v=d s=m*y; w2=1/sqrt(s(3,1) subplot(3,1,2) plot(0 2,0 0) hold on plot(0 0.5 1 1.5 2,0 v(1,1) v(2,1) v(3,1) 0) q=v*v*m; r=j-(s(3,1)/p)*q; d=1 -1
25、1; % x(1)=d; for i=1:16; y=r*d; % l=y(i) d=(1/y(3,1)*y; end,第7章 振动,v=d s=m*y; w3=1/sqrt(s(3,1) subplot(3,1,3) plot(0 2,0 0) hold on plot(0 0.5 1 1.5 2,0 v(1,1) v(2,1) v(3,1) 0),运行结果,第7章 振动,四、 强迫振动及减振器,第7章 振动,m=50 10; k=200 40; c=10 6; N=m(2) c(2) k(2);c(2) k(2); D=m(1)*m(2) (c(1)+c(2)*m(2)+c(2)*m(1)
26、 .(k(1)+k(2)*m(2)+k(2)*m(1)+c(1)*c(2) .k(1)*c(2)+c(1)*k(2) k(1)*k(2); sy=tf(N,D); % y,t=impulse(transferab(m,k,c),20) y,t=impulse(sy,20); subplot(2,1,1); plot(t,y(:,1,1); ylabel(x_1(t); title(impulse response of m_1) subplot(2,1,2); plot(t,y(:,2,1); xlabel(time t); ylabel(x_2(t) title(impulse respon
27、se of m_2);,编制MATLAB文件,运 行 结 果,第7章 振动,减振器：求m1上施加一外力时，物体m1和m2位移的频率响应函数。程序如上,k=200 40; c=10 6; omega=0.0:0.005:4 for i=1:2if i=1sys=tf(1,m(1) c(1) k(1);mag,phas=bode(sys,omega);plot(omega,mag(1,:),-);hold on;else N=m(2) c(2) k(2);c(2) k(2); D=m(1)*m(2) (c(1)+c(2)*m(2)+c(2)*m(1) .(k(1)+k(2)*m(2)+k(2)*m
28、(1)+c(1)*c(2) .k(1)*c(2)+c(1)*k(2) k(1)*k(2); sys=tf(N,D);,sys=tf(N,D);mag,phas=bode(sys,omega);plot(omega,mag(1,:); end end xlabel(excitation frequency(rad/s) ylabel(|x_1|) text(2.1,0.045,no vibration absorber) text(0.3,0.02,with vibration absorber),第7章 振动,运行结果,第7章 振动,习题1 如图示系统求对,的响应,解：,图 7-1,第7章 振
29、动,本章习题,下面的习题已给出了应用振动理论求解的步骤，请同学们参照其求解步骤，用MATLAB编程求解,（1）,初始条件,第7章 振动,振动微分方程,固有频率和主振型,解得,第7章 振动,令,振型矩阵,主刚度矩阵,主质量矩阵,第7章 振动,正则振型矩阵,令 ，方程（1）左乘 后得,第7章 振动,第7章 振动,第7章 振动,习题2 写出用矩阵迭代法求正定系统的固有频率，主振型的步骤并举一两自由度系数进行迭代计算。,解：1.求出系统的质量矩阵M，刚度矩阵K（或柔度矩阵F），计算出动力矩阵A=FM=K-1M,2.任选一非零初始向量X1，最后一个元素为13.作矩阵迭代Y=AX4.重复迭代步骤Y=AX5.若X与X在误差范围内，则迭代结束,第7章 振动,6.迭代r+1阶，用A代替A，重前步骤即可,例如：,第7章 振动,解：,质量矩阵,刚度矩阵,第7章 振动,设,第7章 振动,第7章 振动,第7章 振动,第7章 振动,迭代第二阶,假设初始向量,第7章 振动,本章结束！,
