单盘转子故障诊断报告(附MATLAB程序)

使用Matlab技术进行故障检测的基本步骤随着科学技术的不断发展，人们对于工程设备使用中的故障检测越来越重视。

而Matlab作为一种强大的数学建模和仿真软件，被广泛应用于各个领域中，包括故障检测。

本文将介绍使用Matlab技术进行故障检测的基本步骤，旨在帮助读者了解如何利用该软件来提高设备的可靠性和性能。

1. 数据采集与预处理故障检测的第一步是数据采集。

通过传感器、仪表或其他相关设备，收集机器或系统运行过程中所产生的数据。

这些数据包括振动信号、电流信号、温度信号等，以及与工程设备相关的操作参数。

数据质量的好坏直接影响到后续故障检测的准确性和可信度。

在数据采集之后，对数据进行预处理是必不可少的。

这包括数据滤波、降噪、归一化等过程。

例如，通过将信号进行低通滤波，可以降低高频噪声对故障检测结果的干扰。

通过归一化处理，可以消除由于不同传感器测量范围不同而带来的问题。

预处理的目的是提高原始数据的质量和可用性，为后续的故障检测算法提供更好的输入。

2. 特征提取与选择在完成数据预处理之后，下一步是从原始数据中提取有用的特征。

特征是指反映设备状态和特性的数值或统计量。

常用的特征包括振幅谱、频率谱、时频分析等。

它们能够揭示设备在不同时间和频率上的特定行为，如共振、共振频率变化等。

特征的选择要根据具体的故障类型和设备特性来确定，以充分反映设备的工作状态。

值得一提的是，特征工程是整个故障检测流程中最重要的环节之一。

一个恰当选择的特征集合能够提高故障检测的准确性和鲁棒性。

而Matlab提供了丰富的工具包和函数，以帮助用户进行特征提取和选择工作。

读者可以根据具体需求选择合适的方法和工具来完成特征工程。

3. 故障检测算法的建模与训练在得到有效的特征之后，下一步是选择合适的故障检测算法。

常见的算法有支持向量机（SVM）、人工神经网络（ANN）、卷积神经网络（CNN）等。

这些算法能够根据输入的特征和标签进行模型的建立和训练，从而实现故障检测的功能。

故障诊断第一次实验报告——转子实验台振动信号的时域、频域分析1.由时域信号得出震动的平均值、方差、峭度信息程序代码如下：y=Data(1:6000,1);y=y';[a,b]=xcorr(y,'unbiased');figure(1);plot(b,a);grid;xlabel('位移信号自相关')figure(2);x=0:1/2560:1;plot(x,mean(y));gtext('平均值');hold on;figure(3);plot(x,var(y));gtext('方差');hold on;figure(4);plot(x,kurtosis(y));gtext('峭度');位移信号当转速为600r/min时当转速为1000r/min时当转速为1500r/min时速度信号速度为600速度1000速度15001.由频域信号得出信号的傅里叶变换、功率谱密度信息程序代码如下：t=0:1/2560:1;f=Data;f=f';y=fft(f,82944);m=abs(y);f1=(0:length(y)/2-1)'*2560/length(y); %计算变换后不同点对应的幅值figure(1);plot(f1,m(1:length(y)/2));ylabel('幅值的模');xlabel('时间(s)');title('原始信号傅里叶变换');grid;p=y.*conj(y)/82944; %计算功率谱密度ff=(0:length(y)/2-1)'*2560/length(y); %计算变换后不同点对应的频率值figure(2);plot(ff,p(1:length(y)/2));ylabel('幅值');xlabel('频率(Hz)');title('功率谱密度');grid;位移信号速度600速度1000速度1500速度信号速度600速度1000速度1500。

matlab 模拟电机故障-回复模拟电机故障是现代工程中常见的问题。

正常工作的电机可通过使用MATLAB进行建模和仿真来检查其性能。

然而，当电机发生故障时，需要进行模拟以了解故障的原因和影响。

在本文中，我将详细介绍如何利用MATLAB模拟电机故障。

第一步是建立电机模型。

电机模型可以包括转子、定子、磁路和绕组等组成部分。

可以使用MATLAB中提供的电气特性和动力学建模工具箱来快速建立电机模型。

对于转子部分，可以使用d轴和q轴的电磁方程建立模型。

对于定子部分，可以使用麦克斯韦方程组建立模型。

磁路可以通过磁路方程和调制函数进行建模。

绕组可以通过差分方程和运动方程进行建模。

第二步是定义故障类型。

电机故障可以分为不同的类型，如断相故障、绝缘故障、短路故障等。

每种故障类型会对电机的性能产生不同的影响。

在MATLAB中，可以通过改变电机模型中的特定参数来模拟不同的故障类型。

例如，对于断相故障，可以将相参数设置为零。

对于绝缘故障，可以改变绝缘层的电阻和电容值。

对于短路故障，可以改变绕组的导体尺寸和材料。

第三步是分析故障对电机性能的影响。

使用MATLAB进行模拟后，可以得到电机在不同故障条件下的性能曲线。

可以分析电机的输出功率、效率、转矩和转速等参数的变化情况。

还可以分析电机的失效频率和失效机制。

通过这些分析，可以了解故障对电机的影响程度和持续时间。

第四步是评估故障检测和诊断方法。

在现实世界中，及时检测和诊断电机故障非常重要。

使用MATLAB模拟电机故障可以评估不同的故障检测和诊断方法的有效性。

可以比较不同方法的检测准确性、检出率和误报率等指标。

还可以评估不同故障类型的故障诊断方法的适用性。

第五步是提出故障修复方案。

通过对电机故障的模拟和分析，可以得出修复故障所需的步骤和措施。

例如，对于断相故障，可以建议更换损坏的线圈或重新连接电路。

对于绝缘故障，可以建议加强绝缘层或更换绝缘材料。

对于短路故障，可以建议修复或替换损坏的绕组。

A.1 传递距阵法分析程序%main_critical.m%该程序使用Riccati传递距阵法计算转子系统的临界转速及振型%本函数中均采用国际单位制% 第一步：设置初始条件（调用函数shaft_parameters）%初始值设置包括：轴段数N，搜索次数M%输入轴段参数：内径d,外径D，轴段长度l，支撑刚度K，单元质量mm，极转动惯量Jpp[N,M,d,D,l,K,mm,Jpp]=shaft_parameters;% 第二步：计算单元的5个特征值（调用函数shaft_pra_cal）%单元的5个特征值：%m_k：:质量%Jp_k：极转动惯量%Jd_k：直径转动惯量%EI：弹性模量与截面对中性轴的惯性矩的乘积%rr：剪切影响系数[m_k,Jp_k,EI,rr]=shaft_pra_cal(N,D,d,l,Jpp,mm);% 第三步：计算剩余量（调用函数surplus_calculate）,并绘制剩余量图%剩余量：D1for i=1:1:Mptx(i)=0;pty(i)=0;endfor ii=1:1:Mwi=ii/1*2+50;[D1,SS,Sn]=surplus_calculate(N,wi,K,m_k,Jp_k,JD_k,l,EI,rr);D1;pty(ii)=D1;ptx(ii)=w1endylabel(‘剩余量’)；plot(ptx,pty)xlabel(‘角速度red/s’);grid on% 第四步：用二分法求固有频率及振型图%固有频率：Critical_speedwi=50;for i=1:1:4order=i[D1,SS,Sn]=surplus_calculate(N,wi,k,m_k,Jp_k,Jd_k,l,EI,rr);Step=1;D2=D1;kkk=1;while kkk<5000if D2*D1>0wi=wi+step;D2=D1;[D1,SS,Sn]=surplus_calculate(N,wi,K,m_k,Jp_k,Jd_k,l,EI,rr);endif D1*D2<0wi=wi-step;step=step/2;wi=wi+step;[D1,SS,Sn] =surplus_calculate(N,wi,K,m_k,Jp_k,Jd_k,l,EI,rr);EndD1;Wi;If atep<1/2000Kkk=5000;endendCritical_speed=wi/2/pi*60figure;plot_mode(N,l,SS,Sn)wi=wi+2;end%surplus_calculate,.m%计算剩余量%（1）计算传递矩阵%（2）计算剩余量function [D1,SS,Sn]= surplus_calculate(N,wi,K,m_k,Jp_k,Jd_k,l,EI,rr); % (1)计算传递矩阵%＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝%(a)初值设为0%＝＝＝＝＝＝＝＝＝＝＝＝＝＝＝for i=1:1:N+1for j=1:1:2for k=1:1:2ud11(j,k.i)=0;ud12(j,k.i)=0;ud21(j,k.i)=0;ud22(j,k.i)=0;endendendfor i=1:1:Nfor j=1:1:2for k=1:1:2us11(j,k.i)=0;us12(j,k.i)=0;us21(j,k.i)=0;us22(j,k.i)=0;endendendfor i=1:1:Nfor j=1:1:2for k=1:1:2u11(j,k.i)=0;u12(j,k.i)=0;u21(j,k.i)=0;u22(j,k.i)=0;endendend%＝＝＝＝＝＝＝＝＝＝＝＝%(b)计算质点上传递矩阵―――点矩阵的一部分！%＝＝＝＝＝＝＝＝＝＝＝＝for i=1:1:N+1ud11(1,1,i)=1; ud11(1,2,i)=0; ud11 (2,1,i)=0; ud11(2,2,i)=1;ud21(1,1,i)=0; ud21(1,2,i)=0; ud21 (2,1,i)=0; ud21(2,2,i)=0;ud22(1,1,i)=1; ud22(1,2,i)=0; ud22 (2,1,i)=0; ud22(2,2,i)=1;end%＝＝＝＝＝＝＝＝＝＝＝＝%(c)计算质点上传递矩阵―――点矩阵的一部分！%＝＝＝＝＝＝＝＝＝＝＝＝for i=1:1:N+1ud12(1,1,i)=0;ud12(1,2,i)=(Jp_k(i)-Jd_k(i))*wi^2; %%%考虑陀螺力矩ud12(2,1,i)=m_k(i)*wi^2-k(i);ud12(2,2,i)=0;end%＝＝＝＝＝＝＝＝＝＝＝＝%(d)以下计算的是无质量梁上的传递矩阵―――场矩阵%计算的锥轴的us是不对的，是随便令的，在后面计算剩余量时，zhui中会把错误的覆盖掉%＝＝＝＝＝＝＝＝＝＝＝＝for i=1:1:Nus11(1,1,i)=1; us11(1,2,i)=1(i); us11 (2,1,i)=0; us11(2,2,i)=1;us12(1,1,i)=0; us12(1,2,i)=0; us12 (2,1,i)=0; us12(2,2,i)=0;us21(1,1,i)=1(i)^2/(2*EI(i)); us21(1,2,i)=(1(i)^3*(1-rr(i))/(6*EI(i)); us21(2,1,i)=1(i)/EI(i);us21(2,2,i)=1(i)^2/(2*EI(I));us22(1,1,i)=1; us22(1,2,i)=1(i); us22 (2,1,i)=0; us22(2,2,i)=1; end%＝＝＝＝＝＝＝＝＝＝＝＝%此处全为计算中间量%＝＝＝＝＝＝＝＝＝＝＝＝for i=1:1:N＋2Su (1,1,i)=0; Su (1,2,i)=0; Su (2,1,i)=0; Su (2,2,i)=0;Sn(1,1,i)=0; Sn (1,2,i)=0; Sn (2,1,i)=0; Sn(2,2,i)=0;SS (1,1,i)=0; SS (1,2,i)=0; SS (2,1,i)=0; SS (2,2,i)=0;endfor i=1:1:2for j=1:1:2SS1(i,j)=0;Ud11(i,j)=0; Ud12(i,j)=0; Ud21(i,j)=0; Ud22(i,j)=0;Us11(i,j)=0; Us12(i,j)=0; Us21(i,j)=0; Us22(i,j)=0;endend%＝＝＝＝＝＝＝＝＝＝＝＝%(e)调用函数cone_modify修改锥轴的传递矩阵%＝＝＝＝＝＝＝＝＝＝＝＝cone_modify（4，wi）;cone_modify（5，wi）;cone_modify（6，wi）;cone_modify（7，wi）;cone_modify（8，wi）;cone_modify（16，wi）;cone_modify（17，wi）;cone_modify（18，wi）;cone_modify（19，wi）;cone_modify（22，wi）;cone_modify（24，wi）;%＝＝＝＝＝＝＝＝＝＝＝＝%(f)形成最终传递矩阵%＝＝＝＝＝＝＝＝＝＝＝＝%Ud11 Ud12 Ud21 Ud22 为最终参与计算的传递矩阵for i=1:1:Nu11(:,:,i)=us11(:,:,i)*ud11(:,:,i)+us12(:,:,i)*ud21(:,:,i);u12(:,:,i)=us11(:,:,i)*ud12(:,:,i)+us12(:,:,i)*ud22(:,:,i);u21(:,:,i)=us21(:,:,i)*ud11(:,:,i)+us22(:,:,i)*ud21(:,:,i);u22(:,:,i)=us21(:,:,i)*ud12(:,:,i)+us22(:,:,i)*ud22(:,:,i);endu11(:,:,N+1)=ud11(:,:,N+1); u12(:,:,N+1)=ud12(:,:,N+1);u21(:,:,N+1)=ud21(:,:,N+1); u22(:,:,N+1)=ud22(:,:,N+1);for i=1:1:2for j=1:1:2SS1(i,j)=0;endendfor i=1:1:N+1ud11= u11(:,:,i); ud12= u12(:,:,i); ud21= u21(:,:,i); ud22= u22(:,:,i);SS(:,:,:i+1)=( ud11* SS1+ ud12)*inv(ud21* SS1+ ud22);Su(:,:,i)= ud21* SS1+ ud22;Sn(:,:,i)= inv(ud21* SS1+ ud22); %计算振型时用到SS1=SS(:,:,i+1);end%＝＝＝＝＝＝(2)计算剩余量＝＝＝＝＝＝D1=det(SS(:,:,N+2);for i=1:1:N+1D1=D1*sign(det(Su(:,:,i)); %消奇点end%＝＝＝＝＝＝(2)不平衡响应值EE＝＝＝＝＝＝EE（:,:,n+2）=-inv(SS(:,:,N+2)*PP(:,:,N+2);for i=N+1:-1:1EE(:,:,I)=Sn(:,:,i)*EE(:,:,i+1)-Sn(:,:,i)*UF(:,:,i);endA.2 碰摩转子系统计算仿真程序%main.m%该程序主要完成完成jeffcott转子圆周碰摩故障仿真%＝＝＝＝＝＝＝＝＝＝＝第一步：设置初始条件%rub_sign:碰摩标志，若rub_sign=0,说明系统无碰摩故障；否则rub_sign=1 %loca: 不平衡质量的位置%loc_rub: 碰摩位置%Famp: 不平衡质量的大小单位为：[g]%wi: 转速单位为：[rad]%r: 偏心半径单位为：[mm]%Fampl: 离心力的大小单位为：[kg,m]%fai: 不平衡量的初始相位[rad]clcclear[rub_sign loca loc_rub Famp wi r Famp1 fai]=initial_conditions% 第二步：设置转子系统的参数值%N：划分的轴段数%density: 轴的密度单位为:[kg/m^3%Ef: 轴的弹性模量单位为：[Pa]%L: 每个轴段的长度单位为：[m]%R: 每个轴段的外半径单位为：[m]%ro: 每个轴段的内半径单位为：[m]%miu: 每个轴段的单元质量[kg/m][N density Ef L R Ro miu]=rotor_parameters% 第三步：设置移动单元质量距阵，移动单元质量距阵，刚度单元质量距阵和陀螺力距距阵%Mst: 移动单元质量距阵%Msr: 移动单元质量距阵%Ks: 刚度单元距阵%Ge: 陀螺力距单元距阵[Mst: Msr Ks Ge]=Mst_Msr_Ks_Ge(N,density,R,Ro.L.Ef,miu)% 第四步：距阵组集%M：总的质量距阵%K：总的刚度距阵%G: 总的陀螺力矩距阵[M G K]=M_G_K(N,Ef,R,Ro.Mst,Msr,Ge,Ks,miu,L)% 第五步：加入支撑刚度和阻尼[K C D Ax]=K_D(N,K,M,G)% 第六步：用Newmark方法进行计算%Fen: 每个周期内的步数%ht: 每步的长度ut1=[];xt1=[];yt1=[];N3=(N+1)*4;Fen=100;ht=2*pi/wi/Fen;gama=0.5; beita=0.25;a0=1.0/(beita*ht);al=gama/(beita*ht);a2=1.0/(beita*ht);a3=0.5/beita-1.0;a4=gama/beita-1.0;a5=ht/2.0*(gama/beita-2.0);a6=ht*(1.0-gama);a7=gama*ht;for i=1:1:N3F(i,1)=0;endfor i=1:1;N3yn(i:1)=0;dyn(i:1)=0;ddyn(i:1)=0;endt=0;for n=1:1:Fen*80t=t+ht;n;for i=1:1:30newmark_newton_multiendif mod(n,100)==1nendif n>Fen*60for i=1:1:N+1x(i,1)=yn(i*4-3,n)*le6;y(i,1)=yn(i*4-2,n)*le6;sitax(I,1)=yn(i+4-0,n)/pi*180endu=F(loca*4-4+1,1);ut1=[ut1;[t u]];xt1=[xt1;[t x’]];yt1=[yt1;[t y’]];endendrub_signsave ’xt1.dat’ xt1 –ascii;save ’yt1.dat’ yt1 –ascii;save ’ut1.dat’ ut1 –ascii;%initial_conditions.m%该程序主要进行仿真条件设置Function [rub_sign loca loc_rob Famp wi r Famp1 fai]=initial_conditions%需要设置的初始条件有：%rub_sign: 碰摩标志，若rub_sign=0,说明系统无碰摩故障；否则rub_sign=1 %loca: 不平衡质量的位置%loc_rub: 碰摩位置%Famp: 不平衡质量的大小单位为：[g]%wi: 转速单位为：[rad]% r：偏心半径单位为:[mm]%Famp1：离心力的大小单位为:[kg,m]%fai: 不平衡量的初始相位[rad]rub_sign=0;loca=6;loc_rub=8;Famp=[1];wi=3000/60*2*pi;r=30Famp1=Famp(1)/1000*wi^2*r/1000;fai=[30 30]/180*pi%roto_parameters.m%该程序对Jeffcott转子系统进行参数设置function [N density Ef L R R0 miu]=rotor_parameters%N: 划分的轴段数%density: 轴的密度单位为：[kg/m^3]%Ef: 轴的弹性模量单位为:[Pa]%L 每个轴段的长度单位为:[m]%R 每个轴段的外半径单位为:[m]%R0: 每个轴段的内半径单位为:[m]%miu: 每个轴段的单元质量[kg/m]N=11;Density=7850;Ef=[2.1 2.1 2.1 2.1 2.1 2.1 2.1 2.1 2.1 2.1 2.1]*lell;L=[90.5 90.5 80.5 62.5 30 20 22.5 62.5 90.5 90.5 90.5]/1000;R=[20 20 20 20 20 90 20 20 20 20 20]/2000;R0=[0 0 0 0 0 0 0 0 0 0 0]/2000;for i=1:1;Nmiu(i)=density*pi*(R(i)^2-R0(i)^2)end%Mst_Msr_Ks_Ge.m%该程序设置单元距阵Function [Mst Msr Ks Ge]= Mst_Msr_Ks_Ge(N,density,R,R0,L,Ef,miu)%Mst: 移动单元质量距阵%Msr: 转动单元质量距阵%Ks: 刚度单元距阵%Ge: 陀螺力距单元距阵NN0=N;NN1=NN0+1NN2=NN1+1for i=1:1:NN0Mst(1,1,i)=156;Mst(2,1,i)=0; Mst(2,2,i)=156;Mst(3,1,i)=0; Mst(3,2,i)=-22*L(i); Mst(3,3,i)= 4*L(i)^2;Mst(4,1,i)22*L(i); Mst(4,2,i)=0; Mst(4,3,i)=0;Mst(4,4,i)=4*L(i)^2; Mst(5,1,i)=54; Mst(5,2,i)=0Mst(5,3,i)=0; Mst(5,4,i)=13*L(i); Mst(6,1,i)=0;Mst(6,2,i)=54; Mst(6,3,i)=13*L(i); Mst(6,4,i)=0;Mst(7,1,i)=0; Mst(7,2,i)=13*L(i); Mst(7,3,i)=-3*L(i)^2;Mst(7.4.i)=0; Mst(8,1,i)=-13*L(i); Mst(8,2,i)=0;Mst(8,3,i)=0; Mst(8,4,i)=-3*L(i)^2; Mst(5,5,i)=156;Mst(6,5.i)=0; Mst(6,6,i)=156;Mst(7,5,i)=0; Mst(7,6,i)=22*L(i); Mst(7,7,i)=4*L(i)^2;Mst(8,5,i)=-22*L(i); Mst(8,6,i)=0 Mst(8,7,i)=0Mst(8,8,i)=4*L(i)^2;endfor i=1:1:NN0Msr(1,1,i)=36;Msr(2,1,i)=0; Msr(2,2,i)=36;Msr(3,1,i)=0 Msr(3,2,i)=-3*L(i); Msr(3,3,i)=4*L(i)^2;Msr(4,1,i)=3*L(i); Msr(4,2,i)=0; Msr(4,3,i)=0;Msr(4,4,i)=4*L(i)^2; Msr(5,1,i)=36; Msr(5,2,i)=0;Msr(5,3,i)=0; Msr(5,4,i)=-3*L(i); Msr(6,1,i)=0;Msr(6,2,i)=-36; Msr(6,3,i)= 3*L(i); Msr(6,4,i)=0;Msr(7,1,i)=0; Msr(7,2,i)=-3*L(i); Msr(7,3,i)=-L(i)^2;Msr(7,4,i)=0; Msr(8,1,i)=3*L(i); Msr(8,2,i)=0;Msr(8,3,i)=0; Msr(8,4,i)=-L(i)^2; Msr(5,5,i)=36;Msr(6,5,i)=0; Msr(6,6,i)=36; Msr(7,5,i)=0;Msr(7,6,i)=3*L(i); Msr(7,7,i)=4*L(i)^2; Msr(8,5,i)=-3*L(i);Msr(8,6,i)=0; Msr(8,7,i)=0; Msr(8,8,i)=4*L(i)^2; endfor i=1:1:NN0Ge(1,1,i)=0;Ge(2,1,i)=36; Ge(2,2,i)=0;Ge(3,1,i)=-3*l(i); Ge(3,2,i)=0; Ge(3,3,i)=0;Ge(4,1,i)=0; Ge(4,2,i)=-3*L(i); Ge(4,3,i)=4*L(i)^2;Ge(4,4,i)=0; Ge(5,1,i)=0; Ge(5,2,i)=36;Ge(5,3,i)=-3*L(i); Ge(5,4,i)=0; Ge(6,1,i)=-36;Ge(6,2,i)=0; Ge(6,3,i)=0; Ge(6,4,i)=-3*L(i0);Ge(7,1,i)=-3*L(i); Ge(7,2,i)=0; Ge(7,3,i)=0;Ge(7,4,i)=L(i)62; Ge(8,1,i)=0; Ge(8,2,i)=-3*L(i);Ge(8,3,i)=-L(i)^2; Ge(8,4,i)=0; Ge(5,5,i)=0;Ge(6,5,i)=36; Ge(6,6,i)=0; Ge(7,5,i)=3*L(i);Ge(7,6,i)=0; Ge(7,7,i)=0; Ge(8,5,i)=0;Ge(8,6,i)=3*L(i); Ge(8,7,i)=4*L(i)^2; Ge(8,8,i)=0;endfor i=1:1:NN0Ks(1,1,i)=12;Ks(2,1,i)=0; Ks(2,2,i)=12;Ks(3,1,i)=0; Ks(3,2,i)=-6*L(i); Ks(3,3,i)=4*L(i)^2;Ks(4,1,i)=6*L(i); Ks(4,2,i)=0; Ks(4,3,i)=0;Ks(4,4,i)=4*L(i)^2; Ks(5,1,i)=-12; Ks(5,2,i)=0;Ks(5,3,i)=0; Ks(5,4,i)=-6*L(i); Ks(6,1,i)=0;Ks(6,2,i)=-12; Ks(6,3,i)=6*L(i); Ks(6.4.i)=0;Ks(7,1,i)=0; Ks(7,2,i)=-6*L(i); Ks(7,3,i)=2*L(i)^2;Ks(7,4,i)=0; Ks(8.1,i)=6*L(i); Ks(8,2,i)=0;Ks(8,3,i)=0; Ks(8,4,i)=2*L(i)62; Ks(5,5,i)=12;Ks(6,5,i)=0; Ks(6,6,i)=12; Ks(7,5,i)=0;Ks(7,6,i)=6*L(i); Ks(7,7,i)=4*L(i)^2; Ks(8,5,i)=-6*L(i);Ks(8,6,i)=0; Ks(8,7,i)=0; Ks(8,8,i)=4*L(i)^2; endfor i=1:1:NN0for j=1:1:8for k=1:1:8EI=Ef(i)*pi*(R(i)^4-R0(i)^4-)/4;Mst(j,k,i)=Mst(j,k,i)*miu(i)*L(i)/420;Msr(j,k,i)=Msr(j,k,i)*miu(i)*R(i)^2/120/L(i);Ge(j,k,i)=-Ger(j,k,i)*2*miu(i)*R(i)^2/120/L(i);Ks(j,k,i)=Ks(j,k,i)*EI/L(i)^3;endendendfor i=1:1:NN0for j=1:1:8for k=j:1:8Mst(j,k,i)=Mst(k,j,i);Msr(j,k,i)=Msr(k,j,i);Ks(j,k,i)=Ks(k,j,i);Ge(j,k,i)=-Ge(k,j,i);endendend%M_G_K.m%该程序进行距阵组集function [M G K]=M_G_K(N,Ef,R,R0,Mst,Msr,Ge,Ks,miu,L)% M: 总的质量距阵% K：总的刚度距阵% G：总的陀螺力距距阵NN0=Nfor i=1:1:NN0for j=1:1:8for K=1:1:8Ms(j,k,i)=Mst(j,k,i)+Msr(j,k,i);Ks(j,k,i)=Ks(j,k,i);Ge(j,k,i)=-Ge(j,k,i);endendengfor i=1:1:Nfor j=1:1:8for K=1:1:8M(i*4+j-4,i*4+k-4)=Ms(j,k,i);G(i*4+j-4,i*4+k-4)=Ge(j,k,i);K(i*4+j-4,i*4+k-4)= K (j,k,i);endendendfor i=2:1:Nfor j=1:1;4for k=1:1:4M(i*4+j-4,i*4+k-4)= M(i*4+j-4,i*4+k-4)+Ms(j+4,k+4,i-1);G(i*4+j-4,i*4+k-4)= G(i*4+j-4,i*4+k-4)+Ge(j+4,k+4,i-1);K(i*4+j-4,i*4+k-4)= K(i*4+j-4,i*4+k-4)+Ks(j+4,k+4,i-1);endendendKzx(1)=7e7;Kzy(1)=7e7;Kzx(12)=7e7;Kzy(12)=7e7;for i=1:1:N+1K(i*4+1-4,i*4+1-4)= K(i*4+1-4,i*4+1-4)+Kzx(i);K(i*4+2-4,i*4+2-4)= K(i*4+2-4,i*4+2-4)+Kzy(i);end%K_D.m%该程序将组尼加入总体距阵function [K C D Ax]=K_D(N,K,M,G)Kzx(1)=7e7;Kzy(1)=7e7;Kzx(12)=7e7;Kzy(12)=7e7;for i=1:1:N+1K(i*4+1-4,i*4+1-4)= K(i*4+1-4,i*4+1-4)+Kzx(i);K(i*4+2-4,i*4+2-4)= K(i*4+2-4,i*4+2-4)+Kzy(i);endformat long g[Ax WW]=eig(inv(M)*K);f=sqrt(WW)/(2*pi);f0=diag(f);f00=abs(sort(f0))fn123=[f00(1) f00(3) f00(5)]wi1=54*2*pi;wi2=284*2*pi;yita1=0.1;yita2=0.2;alf=2*(yita2/wi2-yita1/wi1)*(1/wi1^2-1/wi1^2);beita=2*(yita2*wi2-yita1*wi1)/(1/wi1^2-wi1^2);C=alf*M+beita*K;D=C+G;Dzx(1)=7e3;Dzy(1)=7e3;Dzx(12)=7e3;Dzy(12)=7e3;for i=1:1:N+1D(i*4+1-4,i*4+1-4)=D(i*4+1-4,i*4+1-4)+Dzx(i);D(i*4+2-4,i*4+2-4)=D(i*4+2-4,i*4+2-4)+Dzy(i);endD=D*1;%newmark_newton_multi.m%该程序为newmark_newton算法xn=yn(:,n);dxn=dyn(:,n);ddxn=ddyn(:,n);if i==1xn=yn(:,n);endif i>1xn=yn(:,n+1);endK0=K;K1=K*0;F(loca*4-4+1,1)=Famp1*cos(wi*t+fai(1));F(loca*4-4+2,1)=Famp1*sin(wi*t+fai(1));F(loca*4+1,1)=Famp1*cos(wi*t+fai(1));F(loca*4+2,1)=Famp1*sin(wi*t+fai(1));Pn1=F;Fr=Pn1*0;dert=20/1e6;k_rub=5e4;miu_rub=0.2;xx=yn(loc_rub*4-4+1,n)+5/1e6;yy=yn(loc_rub*4-4+2,n);rr=sqrt(xx^2+yy^2);if n>Fen*60&rr>dertrub_sign=1;K1(loc_rub*4-4+1,loc_rub*4-4+1)=k_rub*(rr-dert);K1(loc_rub*4-4+2,loc_rub*4-4+1)=k_rub*(rr-dert)*miu_rub;K1(loc_rub*4-4+1,loc_rub*4-4+2)=-k_rub*(rr-dert)*miu_rub;K1(loc_rub*4-4+2,loc_rub*4-4+2)=k_rub*(rr-dert);Fr(loc_rub*4-4+1,l)=k_rub*(rr-dert)/rr*(xx-miu_rub*yy);Fr(loc_rub*4-4+2,l)=k_rub*(rr-dert)/rr*(yy+miu_rub*xx); endKT=K0+K1;If i==1Fn1=K0*xn+Fr;endif i>1Fn1=K0*xn1+Fr;endKj=KT+a0*M+a1*C;Pj=Pn1+M*(a0*xn+a2*dxn+a3*ddxn)+C*(a1*xn+a4*dxn+a5*ddxn); If i==1Pj=Pj-Fn1+KT*xn;endif i>1Pj=Pj-Fn1+KT*xn1;endif rr<derti=100;end[QQ,RR]=qr(Kj);xn1=RR\QQ’ *Pj;ddxn1=a0*(xn1-xn)-a2*dxn-a3*ddxn;dxn1=dxn+a6*ddxn+a7*ddxn1;yn(:,n+1)=xn1;dyn(:,n+1)=dxn1;ddyn(;,n+1)=ddxn1;A.3 稳定性分析程序% 第一步：设置初始条件[N Fen y d y2]=initial_conditionsfm=[];for kk=1:1:60w=40+(kk-1)*1/1 %频率区间h=1/w/Fen; %每个周期内积分点数为Fen=200 w=2*pi*w;% 第二步：通过打靶法计算振幅的分岔特性shooting;% 第三步：计算floquet稳定性并保存计算数据floquet_stabilityend%shooting.m%该程序主要通过打靶法计算振幅的分岔特征并存储计算数据for i=1:1:30*Fent=0;t=t+(i)*h;y=rkutta(t,h,y,w);if i>(20*Fen) & mod(i,Fen)==1fprintf(f1,w/2/pi,y(1),y(2),y(3),y(4), y(5),y(6),y(7),y(8));endend%floquet_stability.m%该程序主要计算floquet稳定性并保存计算数据y2=eye(N);s0=y;for i=1:1:Fen*1t=0;t=t+(i)*h;y=rkutta(t,h,y,w);At=fun_At(t,y,w); % 根据当前y值求AFor j=1:1:N % 由dy2=A*y2积分求得y2zzzz=rkutta21(t,h,y2(:,j)’,At);endends=y;Rsi=s0-s;S=y2;Drds=eye(N)*1-S;dsi=Drds\(-Rsi)’ ;y=s-dsi’ ;floquet_mul=max(abs(eig(y2)));fm=[fm;w/2/pi floquet_mul];fprintf(f2,’ %15.9f, %15.9f\n’ ,w/2/pi, floquet_mul);% rkutta.mfunction yout=rkutta(t,h,y,w)N=length(y);for i=1:1:Na(i)=0;d(i)=0;b(i)=0;y0(i)=0;enda(1)=h/2; a(2)=h/2;a(3)=h; a(4)=h;d=fun(t,y,w);b=y;y0=y;for k=1:1:3for i=1:1:Ny(i)=y0(i)+a(k)*d(i);b(i)=b(i)+a(k+1)*d(i)/3;endtt=t+a(k);d=fun(tt,y,w);endfor i=1:1:Nyout(i)=b(i)+h*d(i)/6;end%rkutta21.mfunction yout=ekutta21(t,h,y,A)N=length(y);For i=1:1:Na(i)=0;d(i)=0;b(i)=0;y0(i)=0;enda(1)=h/2; a(2)=h/2;a(3)=h; a(4)=h;d=fun21(t,y,A);b=y;y0=y;for k=1:1:3for i=1:1:Ny(i)=y0(i)+a(k)*d(i);b(i)=b(i)+a(k+1)*d(i)/3;endtt=t+a(k);d=fun21(t,y,A);endfor i=1:1:Nyout(i)=b(i)+h*d(i)/6;end%fun.mfunction d=fun(t,y,w)N=length(y); %圆盘1的质量m1=0.7902; %圆盘2的质量m2=0.7902; %转轴的直径d0=0.01;l1=0.16;l2=0.16;l3=0.16;e1=5.695e-3*0; %圆盘1的偏心e2=5.695e-5; %圆盘2的偏心E=2.11e11; %转轴的弹性模量delta1=200e-5; %碰摩间隙1delta2=2e-5; %碰摩间隙2kr1=1e5; %碰摩刚度1kr2=1e5; %碰摩刚度2fr1=0.1; %碰摩摩擦系数1fr2=0.1; %碰摩摩擦系数2I=pi*d0*d0*d0*d0/64;k11=2.3704e5;k22=2.3704e5;k12=-2.0741e5;k21=k12; % 刚度c11=53.118385200000006; c12=-37.333800000000004;c21=-37.333800000000004; c22=53.118385200000006; % //阻尼系数omega=w;g=0;% 转子动力学方程If y(5)>=deltalDirac1=1;endif y(5)<deltalDirac1=0;endif y(7)>delta2Dirac2=1endif y(7)>delta2Dirac2=0;endfx1=-kr1*(y(5)-delta1)*Dirac1;fx2=-kr2*(y(7)-delta2)*Dirac2;fy1=-fr1*kr1*(y(5)-delta1)*Dirac1;fy2=-fr2*kr2*(y(7)-delta2)*Dirac2;d(1)=y(2);d(2)=-(c11/m1)*y(2)-(c12/m1)*y(4)-(k11/m1)*y(1)-(k12/m1)*y(3) +e1*omega*omega*sin(omega*t)+fy1/m1-g;d(3)=y(4);d(4)=-(c21/m2)*y(2)-(c22/m2)*y(4)-(k21/m2)*y(1)-(k22/m2)*y(3) +e2*omega*omega*sin(omega*t)+fy2/m2-g;d(5)=y(6);d(6)=-(c11/m1)*y(6)-(c12/m1)*y(8)-(k11/m1)*y(5)-(k12/m1)*y(7) +e1*omega*omega*cos(omega*t)+fy1/m1;d(7)=y(8);d(8)=-(c21/m2)*y(6)-(c22/m2)*y(8)-(k21/m2)*y(5)-(k22/m2)*y(7+e2*omega*omega*cos(omega*t)+fx2/m2;% fun_At.mfunction At=fun_At(t,y,w)N=length(y);eps=1e-6;for i=1:1:Nyi=y;yi(i)=y(i)+y(i)*eps;At0=fun(t,y,w);Ati=fun(t,yi,w);For j=1:1:NAt(j,i)=(Ati(j)-At0(j))/(yi(i)*eps);endendA.4 不对中转子系统仿真程序%如下是考虑轴承不对中产生的附加载荷。

Matlab中的异常检测和故障诊断方法在工程和科学领域中，异常检测和故障诊断是非常重要的任务。

随着技术的发展和数据量的增加，传统的方法已经无法满足对异常事件和故障的准确识别和定位的需求。

作为一种功能强大且易于使用的计算工具，Matlab为异常检测和故障诊断提供了多种高级方法和技术。

一、异常检测方法在Matlab中，有多种异常检测方法可以选择。

其中一种常用的方法是基于统计学原理的方法，如均值和方差的检测。

这些方法假设数据服从正态分布，通过计算数据的均值和方差与预先设定的阈值进行比较，来判断数据是否异常。

然而，这种方法对于非正态分布的数据效果可能不佳。

另一种常用的异常检测方法是基于机器学习的方法，如聚类、分类和离群点检测。

在Matlab中，可以使用机器学习工具箱来实现这些方法。

例如，可以使用支持向量机（SVM）算法来进行异常检测。

该算法通过训练一个分类器来划分正常和异常数据，然后对新样本进行分类。

如果新样本被分类为异常，则可以认为该样本是异常数据。

除了以上方法，还有一些其他的异常检测方法，如基于神经网络的方法和基于时间序列的方法。

这些方法有各自的特点和适用范围，根据具体问题选择合适的方法进行异常检测。

二、故障诊断方法与异常检测相比，故障诊断更加复杂和困难。

故障诊断是在异常检测的基础上，进一步分析和推断异常的原因和根源。

在Matlab中，也有多种故障诊断方法可以选择。

其中一种常用的故障诊断方法是基于逻辑推理的方法，如专家系统和诊断树。

这些方法通过构建一个基于规则的知识库，根据给定的观测数据和先验知识，推断系统中存在的故障。

例如，可以使用模糊逻辑方法来构建一个模糊规则库，根据输入数据的模糊集合和规则的模糊推理，来诊断系统中的故障。

另一种常用的故障诊断方法是基于统计模型的方法，如隐马尔可夫模型（HMM）和贝叶斯网络。

这些方法通过建立一个概率模型，根据给定的观测数据和系统模型，推断系统中存在的故障。

例如，可以使用HMM来建立一个模型，将系统的状态表示为隐变量，并通过观测数据进行学习和推断，来诊断系统中的故障。

如何在Matlab中进行数据故障诊断数据故障诊断是在数据分析中的一项重要任务。

在现代科技发展的背景下，随着数据的海量积累和快速传输，数据故障诊断显得尤为重要。

其中，Matlab作为一款功能强大的数据处理和分析工具，为我们提供了丰富的数据故障诊断方法和工具。

下面将就如何在Matlab中进行数据故障诊断进行介绍。

首先，我们要了解什么是数据故障诊断。

数据故障诊断是指通过对数据进行分析和处理，识别出其中的异常值、错误或潜在故障，帮助我们发现和解决数据中的问题。

在数据故障诊断中，我们常常需要用到统计方法和机器学习算法来辅助分析。

在Matlab中，我们可以利用一些基本的数据故障诊断方法来进行数据分析。

首先是利用统计方法识别异常值。

Matlab中提供了一系列的统计函数，如mean、std等，可以用来计算数据的均值和标准差，通过这些统计量，我们可以判断数据中是否存在异常值。

此外，Matlab还提供了hist函数用于绘制数据分布直方图，通过观察数据的分布形态，我们也可以大致判断出数据中是否存在异常。

除了统计方法，Matlab还提供了一些经典的机器学习算法，如聚类分析、主成分分析等，用于进行复杂的数据故障诊断。

聚类分析可以将数据集根据数据之间的相似性进行分组，通过观察每个群组内的数据情况，我们可以发现其中的异常群组。

而主成分分析可以将多个变量通过线性组合的方式转化为少数几个主成分，通过分析主成分的贡献度，我们可以了解各个变量对数据的影响程度，从而发现数据中的异常变量。

除了以上提到的数据故障诊断方法，Matlab还提供了一些其他有用的工具和函数。

例如，异常检测工具箱（Anomaly Detection Toolbox）提供了一系列用于异常检测的算法和函数，包括基于统计学的方法、基于机器学习的方法等。

通过使用这些工具和函数，我们可以更加方便地进行数据故障诊断。

此外，Matlab还提供了数据可视化工具，如绘制散点图、曲线图等，可以帮助我们更直观地观察数据的分布和趋势，从而更好地进行故障诊断。

摘要现代机械设备发展的一个明显的趋势是向大型化、高速化、连续化和自动化的方向发展。

由此而导致设备的功能越来越多，性能指标越来越高，组成和结构越来越复杂，同时对设备管理与维修人员的素质要求也越来越高。

保证高效、安全生产，减少设备特别是重要设备因故障引起的停机时间，是企业提高经济效益和社会效益的前提，而设备状态监测与故障诊断技术将为此提供一个有效的解决途径。

在机械学领域内，机械故障诊断技术可谓是一项富有活力的新技术。

近三十年来，随着微电子技术，计算机技术以及由此而引发的信息技术革命向古老传统机械科学领域的渗透，机械故障诊断技术应运而生。

近年来，工业界、学术界对现代机械故障诊断技术认识的深化与日俱增。

本论文以转子为主要研究对象，对其振动监测及故障诊断的技术进行了分析研究。

根据故障诊断的基本工作原理，基于MATLAB系统对转子系统产生的振动信号进行采集、存取并进行时频域分析。

在此基础上更利于旋转机械故障的诊断与处理从而减少机械故障给生产带来的损失关键词：转子；旋转机械；故障诊断；MATLABAbstractThe development of modern machinery and equipment is a cleartrend is the large-scale,high-speed,continuous and automaticdirection.Thus cause the device to more and more functions,more and more high performance,composition and structure of more complex,while the equipment management and maintenance staff,the quality requirements are also ncreasing. Ensure the efficient and safe production, it is important to reduce the equipment, particularly equipment downtime caused by failures, economic and social benefits of enterprise to improve the premise, and equipmentcondition monitoring and fault diagnosis technology would providean effective solution.In the field of mechanics, machine fault diagnosis techniques can be described as a dynamic new technology. Past three decades, with the micro-electronics technology, computer technology, and thustriggered the ancient tradition of the information technology revolution to the penetration of mechanical science, mechanical fault diagnosis technology came into being. In recent years, the industrial sector, academia, fault diagnosis technology of modern machineryand deepens the understanding growing.In this thesis, the rotor as the main research object, its vibration monitoring and fault diagnosis techniques are analyzed. Faultdiagnosis based on the basic working principle of the rotor systembased on MATLAB system vibration signals generated byacquisition, access and the time and frequency domain analysis. On this basis, more conducive to rotating machinery fault diagnosis andtreatment thereby reducing mechanical failure caused the loss to the production.KeyWords: rotor, rotating machinery, fault diagnosis, MATLAB目录摘要 (I)Abstract (I)1 绪论 (1)1.1 课题的研究背景和研究意义 (1)1.2 机械设备诊断技术的意义 (2)1.3 机械设备故障诊断的基本方法 (4)1.4转子诊断技术与应用现状 (5)1.5转子故障诊断技术的发展方向 (7)1.6本课题主要研究内容 (7)本章小结 (8)2 旋转机械典型故障的研究与诊断实例 .............................................................- 9 -2.1 旋转机械的典型故障原因及特征分析 .....................................................- 9 -2.1.1 转子不平衡(失衡) .............................................................................- 9 -2.1.2 转子不对中 .....................................................................................- 12 -2.1.3 基座或装配松动 .............................................................................- 14 -2.1.4 油膜涡动及油膜振荡 .....................................................................- 14 -2.1.5 旋转失速和喘振 .............................................................................- 14 -2.2 电动机转子系统振动监测及故障诊断的技术研究 ...............................- 14 -本章小结 ..........................................................................................................- 16 -3 转子的故障信号分析处理 .. (17)3.1 信号的频域分析方法 (17)3.2 转子故障的判断标准 (21)3.2.1 按轴承振幅的评定标准 (22)3.2.2 按轴承烈度的评定标准 (23)3.2.3 按轴振幅的评定标准 (25)3.2.4 其他有关标准 (26)4 故障诊断系统总体设计 ...................................................................................- 29 -4.1 Matlab软件简介........................................................................................- 29 -4.2 转子故障诊断系统总体设计 ...................................................................- 31 -4.2.1 系统界面子系统 .............................................................................- 31 -4.2.2 数据加载子系统 .............................................................................- 32 -4.2.3 时域分析子系统 .............................................................................- 32 -4.2.4 频域分析子系统 .............................................................................- 33 -4.2.5 打印子系统 .....................................................................................- 35 -5 结论 (36)致谢 (37)参考文献 (38)1 绪论1.1 课题的研究背景和研究意义航空发动机、汽轮机、压缩机、风机、水泵等旋转机械，在国防、能源、电力、交通、机械和化工等领域中广泛应用并发挥着重要作用。

Matlab技术在机器故障诊断中的应用方法在现代制造业中，机器故障诊断一直是一个重要的问题。

当机器发生故障时，能够快速准确地判断故障原因，对于提高生产效率、降低维修成本至关重要。

传统的故障诊断方法通常基于经验和直觉，但这种方法存在一定的主观性和不准确性。

随着计算机技术的不断发展，Matlab技术在机器故障诊断中的应用逐渐得到重视。

Matlab是一种功能强大的数值计算和数据分析软件，它不仅拥有丰富的数学和统计函数库，还提供了强大的绘图和可视化功能。

这使得Matlab成为机器故障诊断的理想工具。

下面将介绍几种基于Matlab的机器故障诊断方法。

首先，基于信号处理的故障诊断方法是最常见的一种。

通过对机器传感器采集到的信号进行处理和分析，可以得到机器的工作状态和故障特征。

例如，通过对机器振动信号的频谱分析，可以检测出机器可能存在的故障，如轴承磨损、不平衡等。

Matlab提供了丰富的信号处理工具箱，可以方便地进行频谱分析、时频分析等操作。

此外，Matlab还支持自定义算法和模型，可以根据实际情况对信号进行特征提取和分类。

其次，基于机器学习的故障诊断方法也是一种有效的手段。

机器学习是一种从数据中自动学习模式和规律的技术，可以应用于故障诊断。

通过对大量样本数据进行训练和学习，可以建立机器故障诊断模型，实现自动识别和分类。

Matlab提供了丰富的机器学习工具箱，如支持向量机、神经网络、决策树等，可以方便地进行模型训练和预测。

此外，Matlab还提供了数据可视化工具，可以对机器学习结果进行直观的展示和分析。

另外，基于图像处理的故障诊断方法也是一种有效的手段。

对于一些机器设备，如电子元件、机器表面等，可以通过图像处理技术来检测和诊断故障。

例如，通过对机器图像进行边缘检测、目标检测等操作，可以得到机器的几何特征和表面损伤情况。

Matlab提供了丰富的图像处理工具箱，可以方便地进行图像增强、边缘检测、目标跟踪等操作。

此外，Matlab还支持图像特征提取和分类，可以根据实际情况构建机器故障诊断模型。

Matlab中的故障诊断与预测方法介绍近年来，随着工业自动化的发展，对设备的故障诊断与预测需求越来越迫切。

Matlab作为一种功能强大的数值计算与数据分析工具，提供了丰富的故障诊断与预测方法。

本文将介绍一些常用的Matlab中的故障诊断与预测方法，包括信号处理、统计分析和机器学习方法。

一、信号处理方法在故障诊断与预测中，信号处理是一种常用的方法。

Matlab提供了丰富的信号处理工具箱，可以实现对信号的滤波、频谱分析和时频分析等操作。

对于故障诊断来说，信号滤波是一项基本任务。

Matlab中可以使用多种滤波器设计方法，如低通、高通、带通和带阻滤波器等。

通过滤波处理，可以去除噪声，提取有用的故障特征。

另外，频谱分析也是一种常用的信号处理方法。

通过对信号进行傅里叶变换，可以将信号从时域转换到频域。

在Matlab中，提供了多种频谱分析函数，如功率谱密度估计、频谱图和频谱分析等。

这些函数可以帮助用户深入分析信号的频谱特征，从而进行故障诊断。

时频分析是信号处理的另一种重要方法。

它可以将信号从时域和频域两个方面进行分析，从而更加全面地了解信号的特征及变化规律。

Matlab中提供了多种时频分析函数，如连续小波变换、峰值分析和瞬时频率分析等。

这些函数可用于分析非平稳信号的时频特性，为故障诊断提供更加准确的判断。

二、统计分析方法除了信号处理方法，统计分析也是故障诊断与预测的重要手段。

Matlab提供了丰富的统计工具箱，包括概率分布、假设检验、方差分析和回归分析等。

概率分布可以用于描述数据的分布情况。

在故障诊断中，通过对实际采集的数据进行概率分布拟合，可以判断故障是否存在。

Matlab中提供了常用的概率分布函数，如正态分布、指数分布和泊松分布等。

这些函数可用于分析数据的分布特征，为故障诊断提供统计依据。

假设检验可以用于验证某个假设是否成立。

在故障诊断中，可以基于假设检验的方法进行故障的判断。

Matlab中提供了多种假设检验函数，如单样本t检验、双样本t检验和方差分析等。