机械振动大作业

1、机械振动项目2017-5-21沈阳航空航天大学学院：创新创业学院班级：班组员1：邵增德学号：1 组员2：哈哲远学号：0目录一、significance of the project2二、Matrix Iteration method22.1 方法的一般性说明或介绍22.2 方法的应用说明32.3 计算结果5三、Transfer matrix method63.1 方法的一般性说明或介绍63.2 方法的应用说明73.3计算结果8四、近似计算结果的验证9五、两种近似方法及其计算结果的比较分析10附录11附录 A The Matlab program of matrix iteration meth

2、od11附录 B The Matlab program of transfer matrix method12一、 significance of the project Aero engine modern pursuit of high performance and high thrust weight ratio, the structure is more complex, more and more harsh working conditions, which leads to the excessive vibration factors gradually increased

3、. Vibration is the decisive factor that affects the life of aero-engine and flight safety. The following figure shows the main shaft of an aero-engine compressor and turbine, which is simplified as a multi degree of freedom vibration system. In order to ensure its safe operation, it is necessary to

4、determine its main vibration modes, natural frequencies and main vibration laws.aaaaaa6aa二、 Matrix Iteration method2.1 方法的一般性说明或介绍系统的主振型方程可分别表示为 （1）。 （2）矩阵迭代法是从假设主振型出发，对以上两式进行矩阵迭代运算。对（1）式进行迭代时，是从系统的最高阶固有频率与主振型开始，逐步求得系统的全部或一部分固有频率和主振型。而对于（2）进行迭代时，则以系统的最低阶固有频率与主振型开始，逐步求得系统的全部或一部分固有频率和主振型。因为工程上通常对系统的最低

5、阶或较低的几阶固有频率和主振型比较重视，所以下面只考虑的（2）式的迭代运算。引进动力矩阵, (3)则（2）变为 (4)当（4）式中的A取为系统的各阶主振型时，该式成为 (5)当（5）式中的A取为系统各阶主振型时，该式成为 （i=1,2.，r） （6）随着迭代次数的增加，AK越来越接近于第一阶主阵型，即,AK=AK-1，因此我们将得到关于固有频率的关系式，从而计算出n阶固有频率。2.2 方法的应用说明1、本方法为矩阵迭代法，所需要的公式如下所示：转动惯量 J=12mR2刚度 K=G*IpL惯性矩 Ip=D4322、根据所给条件，进行编程计算刚度k和转动惯量J，为了计算的方便，程序如下：/ D1=

6、0.4;h=0.04;G=7.69e10;rou=7800;D2=0.04;a=0.12;V=pi*(D12)*0.25*h;m=rou*V;j=0.5*m*(D12)*0.25;Ip=pi*D24*(1/32);k=(G*Ip)/a; /3、对系统各段的转动惯量以及刚度进行赋值并分别列出对应的矩阵，具体程序如下所示：/ J1=j; J2=j; J3=j;J4=j;J5=j;J6=j;J7=j;J8=j;k1=k; k2=k; k3=k;k4=k;k5=k;k6=k;k7=k/6;k8=k;M=J1 0 0 0 0 0 0 0; 0 J2 0 0 0 0 0 0; 0 0 J3 0 0 0 0

7、 0; 0 0 0 J4 0 0 0 0; 0 0 0 0 J5 0 0 0; 0 0 0 0 0 J6 0 0; 0 0 0 0 0 0 J7 0; 0 0 0 0 0 0 0 J8 ;K=k2 -k2 0 0 0 0 0 0; -k2 k2+k3 -k3 0 0 0 0 0; 0 -k3 k3+k4 -k4 0 0 0 0; 0 0 -k4 k4+k5 -k5 0 0 0; 0 0 0 -k5 k6+k5 -k6 0 0; 0 0 0 0 -k6 k6+k7 -k7 0 ; 0 0 0 0 0 -k7 k7+k8 -k8; 0 0 0 0 0 0 -k8 k8 ; /4、然后进行编写矩阵迭

8、代方法的程序，首先根据我们的编程思路画出流程框图如下图2-1所示图2-1具体程序如下所示：/*/ K=K+alf*M;D=inv(K)*M;for j=1:n A=1 0 0 0 0 0 0 0; for i=1:10000 B=D*A; a=B(n); B=B/a; if max(abs(B-A)0. om(j)=1/sqrt(a); om(j)=sqrt(abs(om(j)2-alf); break end A=B; end /*/5、最后计算出结果，并画出主阵型D=D*(eye(n)-B*B*M/(B*M*B); MM(:,j)=B;endMMMMF=0 0 0 0 0 0 0 0; M

9、M; 0 0 0 0 0 0 0 0; plot(MMF)2.3 计算结果1、 根据所编写的程序计算出八阶固有频率如下图2-2所示图2-22、 作出的八阶主阵型图如下图2-3所示：图2-3三、 Transfer matrix method3.1 方法的一般性说明或介绍 对于盘轴扭振系统，可以将该整体系统分割成结构形式相同的单元，即一段轴连接一个圆盘，将盘视为刚性质量元件、轴段视为柔性无质量元件，每个单元的扭转运动取决于轴段梁福安的转角和扭矩，只有这两个主要因素，因此取各元件的转角和扭矩作为状态变量，表示为，对整个盘轴扭振系统可以理解为状态变量是通过各个单元从系统的左侧传递到系统的右侧，我们通过

10、对一个单元进行分析，得出各个单元的传递矩阵，进而可以求出整个系统的传递矩阵，然后再根据边界条件就可以求得关于的方程，由于方程过于复杂，我们可以将该表达式的图形作出来，找出使该表达式为0的点，该点即是我们要求的值，从而求出该系统的各阶固有频率。取单元中的第i个单元来来进行分析：（1）：对刚性质量元件盘进行分析可得：，写成矩阵的形式为：,该矩阵称为点传递矩阵。（2）：对柔性无质量元件轴进行分析可得：，写成矩阵的形式为：， 该矩阵称为场传递矩阵。于是得到第i单元两端状态的传递矩阵（从第i-1个盘右侧到第i个盘右侧），该矩阵为单元传递矩阵。对于有n个盘的轴系，将各个单元的单元传递矩阵反顺序连乘就可得到

11、整个系统两端状态的矩阵形式：,即： 再根据轴系两边的边界条件，在该系统中，可以得到关于的方程，对该关系式作图求出与x轴的交点就可以得到系统的固有频率。3.2 方法的应用说明 1、本方法所使用的的数据和公式同2.2.1步骤，赋值程序也是相同的，故这部分不再陈列出来。2、对系统各段的转动惯量以及刚度进行赋值并分别列出对应的矩阵，具体程序如下所示： J1=j; J2=j; J3=j;J4=j;J5=j;J6=j;J7=j;J8=j;k1=k; k2=k; k3=k;k4=k;k5=k;k6=k;k7=k/6;k8=k;J=J1 J2 J3 J4 J5 J6 J7 J8; k=k1 k2 k3 k4

12、k5 k6 k7 k8;3、然后进行编写传递矩阵方法的程序，首先根据我们的编程思路画出流程框图如下图2-4所示图2-4具体程序如下所示：step=0.1;mn=1;T21p=0;for i=1:9000 om=(i-1)*step; T=eye(2); omx(i)=om; z(i)=0; for j=1:n T=1 1/k(j); -J(j)*om2 1-J(j)*om2/k(j)*T; end T21(i)=T(2,1); if T(2,1)=0 omr(mn)=om; mn=mn+1; end if T(2,1)*T21pn break end T21p=T(2,1);endomr3、

13、最后计算出结果，并画出主阵型 for i=1:n om=omr(i); SV(:,1)=1 0; for j=1:n SV(:,j+1)=1 1/k(j); -J(j)*om2 1-J(j)*om2/k(j)*SV(:,j) end MM(:,i)=SV(1,:);endplot(MM)3.3 计算结果1、根据所编写的程序计算出八阶固有频率如下图2-5所示图2-52、作出的八阶主阵型图如下图2-6所示：图2-6四、 近似计算结果的验证1、 采用特征根特征向量法进行精确计算并验证，具体实现程序如下所示 alf=1;C lam=eig(inv(M)*K);for i=1:n omeig(i)=sq

14、rt(abs(lam(i,i);endfor w=1:n C(:,w)=C(:,w).*(1/C(1,w);endomeigCx=1 2 3 4 5 6 12 13;plot(x,C)grid on 2、 运算结果结果以及图像如下图4.1和图4.2所示图4.1图4.2五、 两种近似方法及其计算结果的比较分析根据传递矩阵和矩阵迭代两种方法的计算，我们不难看出传递矩阵方法计算出来的结果比矩阵迭代法的结果偏大一点，不过几乎是相等的，通过特征值特征向量法的精确计算验证，传递矩阵方法所绘制的主阵型更接近与实际。附录附录 A The Matlab program of matrix iteration m

15、ethod% Tranfer matrixclc; clearn=8;D1=0.4;h=0.04;G=7.69e10;rou=7800;D2=0.04;a=0.12;V=pi*(D12)*0.25*h;m=rou*V;j=0.5*m*(D12)*0.25;Ip=pi*D24*(1/32);k=(G*Ip)/a;J1=j; J2=j; J3=j;J4=j;J5=j;J6=j;J7=j;J8=j;k1=k; k2=k; k3=k;k4=k;k5=k;k6=k;k7=k/6;k8=k;J=J1 J2 J3 J4 J5 J6 J7 J8; k=k1 k2 k3 k4 k5 k6 k7 k8;step=

16、0.1;mn=1;T21p=0;for i=1:9000 om=(i-1)*step; T=eye(2); omx(i)=om; z(i)=0; for j=1:n T=1 1/k(j); -J(j)*om2 1-J(j)*om2/k(j)*T; end T21(i)=T(2,1); if T(2,1)=0 omr(mn)=om; mn=mn+1; end if T(2,1)*T21pn break end T21p=T(2,1);endomrfor i=1:n om=omr(i); SV(:,1)=1 0; for j=1:n SV(:,j+1)=1 1/k(j); -J(j)*om2 1-

17、J(j)*om2/k(j)*SV(:,j) end MM(:,i)=SV(1,:);endplot(MM)附录 B The Matlab program of transfer matrix method% Matrix iteration clc; clearn=8;D1=0.4;h=0.04;G=7.69e10;rou=7800;D2=0.04;a=0.12;V=pi*(D12)*0.25*h;m=rou*V;j=0.5*m*(D12)*0.25;Ip=pi*D24*(1/32);k=(G*Ip)/a;J1=j; J2=j; J3=j;J4=j;J5=j;J6=j;J7=j;J8=j;k1

18、=k; k2=k; k3=k;k4=k;k5=k;k6=k;k7=k/6;k8=k;M=J1 0 0 0 0 0 0 0; 0 J2 0 0 0 0 0 0; 0 0 J3 0 0 0 0 0; 0 0 0 J4 0 0 0 0; 0 0 0 0 J5 0 0 0; 0 0 0 0 0 J6 0 0; 0 0 0 0 0 0 J7 0; 0 0 0 0 0 0 0 J8 ;K=k2 -k2 0 0 0 0 0 0; -k2 k2+k3 -k3 0 0 0 0 0; 0 -k3 k3+k4 -k4 0 0 0 0; 0 0 -k4 k4+k5 -k5 0 0 0; 0 0 0 -k5 k6+k5 -k6 0 0; 0 0 0 0 -k6 k6+k7 -k7 0 ; 0 0 0 0 0 -k7 k7+k8 -k8; 0 0 0 0 0 0 -k8 k8 ;alf=1;C lam=eig(inv(M)*K);for i=1:n omeig(i