一种转子-轴承系统的横向裂纹故障特征量提取方法与流程

本发明属于故障诊断方法
技术领域：
，涉及一种转子-轴承系统的横向裂纹故障特征量提取方法。
背景技术：
转子-轴承系统运行环境复杂，状态多变，有其复杂性和非线性，导致了其在运转中会产生多种类型的故障，有资料表明，裂纹故障占有相当数量的比例，其中又以横向裂纹故障为主。作为旋转机械关键部件的转轴，常在温度多变、压力大、荷载多变甚至有一定腐蚀性的环境中工作，加之制造中常存在的缺陷、运行中摩擦造成损伤、复杂多变的应力场和温度场都易使转轴出现疲劳裂纹，如不能在裂纹产生早期对其进行有效的诊断，裂纹极易迅速扩展，最终导致转轴失稳甚至断裂，造成重大经济损失和人员安全威胁。建立合理的裂纹转子-轴承系统数学模型分析其动力学特性对裂纹故障的识别不但有现实意义，也是裂纹故障诊断的逆问题。自转子动力学成为一门系统的学科以来，裂纹转子的数学建模、动力学特性分析及有关的诊断方法等方面的研究取得了长足的发展。有文献表明横向裂纹具有明显的倍频特征，但实际上转子-轴承通常作为统一的系统存在，而轴承所受非线性油膜力在一定程度上会影响裂纹转子的倍频特征，甚至完全掩盖裂纹转子的特点，即亚临界共振现象消失，常规转子动力学分析方法无法将裂纹转子与正常工作转子区分开来，特别是当裂纹较小时，裂纹的存在无法通过时域图、庞加莱截面图、轴心轨迹图、频谱图体现出来。鉴于此，将裂纹转子-轴承系统整体建模，并引入一种对裂纹变化具有更高敏感性的方法具有现实意义。传统的时频分析方法如傅立叶变换方法，每一时刻的信号都包含了系统的所有频率的振动，每一频率都是全部时间域内所有振动成分的共同作用，因此，其中各种频率成分交杂，相互影响，其对于平稳信号有较好的适应性，但对于非平稳信号无法得到令人满意的效果，而由于裂纹使得系统的刚度随着时间的变化而变化，导致傅立叶分析方法不能很好的描述系统的本质特征。经验模态分析(emd)方法完全跳出了傅立叶方法的限制，其根据数据自身的特性提取信号各个成分，将原始信号分解为从高到低的一系列固有模态函数(imf)，代表了信号内部的本质特征，相比传统方法具有更高的分辨率。与此同时，volterra模型是一种时间序列分析方法，其模型参数凝聚了系统状态的重要信息，准确的volterra模型能够深刻、集中地表达动态系统的客观规律，但是采用volterra模型直接对结构的振动信号进行分析效果并不理想。技术实现要素：本发明的目的是提供一种转子-轴承系统的横向裂纹故障特征量提取方法，将所得系统响应加入随机白噪声模拟实际信号，采用emd-ciit降噪处理后建立不同裂纹角下的emd-volterra模型级数矩阵，得到矩阵的奇异值作为故障特征量，可以分析轴承处在非线性油膜力下不同裂纹角对系统响应的影响。本发明所采用的技术方案是，一种转子-轴承系统的横向裂纹故障特征量提取方法，具体按照以下步骤实施：步骤一：建立裂纹转子-轴承系统动力学模型，获的系统动力学方程，根据系统动力学方程计算不同裂纹角下的系统响应信号；步骤二：对不同裂纹角下的系统响应加入白噪声模拟实际信号，并对加入白噪声后的信号进行emd-ciit降噪处理，得到降噪后的信号；步骤三：对降噪后的信号进行emd分解，提取imf分量；步骤四：对步骤三提取的每个imf分量构造volterra级数模型，得出各imf分量的volterra模型参数后，将各imf分量的volterra模型参数作为矩阵的横矢量构造预测参数矩阵；步骤五：求解所得矩阵的奇异值作为系统的故障特征量。本发明的特征还在于，步骤一中所述的裂纹转子-轴承系统动力学模型包括一转轴以及在转轴左右两端呈对称设置的滑动轴承及其转动设置在转轴中间的质量偏心转子，该裂纹转子-轴承系统所产生的裂纹位于转轴与转子转盘交界处的转轴上。记位于转轴左端滑动轴承的几何中心为o，集中质量为m1，左端滑动轴承在转轴上的中心为o1，x1、y1为左端滑动轴承在转轴x，y方向径向位移，轴承长度l，转子几何中心为o2，集中质量为m2，x2，y2为转子在径向x，y方向位移，e为转子偏心距，滑动轴承的轴套与转轴之间由一层油膜隔开，油膜厚度为h，r为转轴半径，则裂纹转子-轴承系统动力学方程为：其中，k为转子刚度，c1滑动轴承阻尼系数，c2转子阻尼系数，fx、fy分别为滑动轴承在x、y方向所受到的非线性油膜力，g为重力加速度，β为转子偏心相位，ω为转子转速，t为时间。令x2＝x2/c，y2＝y2/c，x1＝x1/c，y1＝y1/c，τ＝(ωt+β)，ε＝e/c，其中，c为滑动轴承的径向间隙，得系统无量纲化方程为：其中，fx和fyfy分别为滑动轴承在x、y方向所受到的无量纲油膜力。转子刚度在固定坐标系o-x-y下的表达式为，式中：式中，k0为完整转子刚度，以裂纹截面圆心o'为原点，垂直横向裂纹的直线为ζ轴，平行横向裂纹的直线为η轴建立旋转坐标系o'-ζ-η，θ为旋转角度，σ为偏心量与旋转坐标轴ζ的夹角，α为裂纹角的一半，δkξ、δkη为ζ、η方向的刚度减少量，为开关函数，无量纲油膜力fx和fy表示为：式中：δ为油膜力的起始角位移，通过进行计算。步骤二的具体过程为：步骤2.1，对不同裂纹角下的系统响应信号加入白噪声模拟实际信号，得到带噪声信号x(t)；步骤2.2，首先确定x(t)所有极值点，用三次样条插值函数分别连接所有极大值、极小值形成上下包络线，x(t)与上下包络线的均值m1的差记为h1，h1＝x(t)-m1h1；步骤2.3，计算x(t)的第一个imf分量c1(t)：若h1满足imf条件，则h1即为x(t)的第一个imf分量c1(t)；若h1不满足imf条件，则将h1视为原始带噪声信号，重复步骤2.2循环q次，直到h1q＝h1(q-1)-m1q满足imf条件，h1q为第q次时的原始带噪声信号，h1(q-1)为第q-1次时的原始带噪声信号，m1q为h1(q-1)与其上下包络线的均值，此时，可得c1(t)＝h1q，c1(t)为x(t)的第一个imf分量；步骤2.4，从x(t)分离出c1(t)，即就是：r1＝x(t)-c1(t)，将r1当成原始带噪声信号重复步骤2.2-步骤2.3得到x(t)的第二个imf分量c2(t)；步骤2.5，将r2当成原始带噪声信号重复步骤2.2-步骤2.3得到x(t)的第三个imf分量c3(t)，其中，r2＝x(t)-c1(t)-c2(t)，依次类推，直到rl为一单调函数不可再分为止，得到l个本征模态函数步骤2.6，通过小波阈值降噪方法将噪声信号和有用信号分离，获得第一层imf分量中的降噪信号并计算c1(t)中的真实噪声信号步骤2.7，对噪声信号c1(n)的采样点进行重新排序，得到新的噪声信号c1(a)(t)和新的原始信号其中，c1(a)(t)＝alter(c1(n))，步骤2.8，对新的原始信号进行emd分解，并对各imf分量采用emd-it降噪处理，得到降噪处理后的各imf分量步骤2.9，获得的各个分量进行emd逆过程，即每个分量相加求和即可得到首个降噪信号x～1(t)；步骤2.10，然后重复步骤2.7-2.9共z-1次，共获得z个降噪信号，分别为步骤2.11，将步骤2.10求得的所有降噪信号进行求均值，得到最终的降噪后的信号imf条件为：(1)在完整信号长度上过零点和极值点的数量相等或不多于一个；(2)任意时刻信号上下包络线关于时间轴对称。步骤四的具体过程为：步骤4.1，计算步骤三提取的各个imf分量的volterra预测系数；步骤4.2，得出各imf分量的volterra模型参数后，将其作为矩阵的横矢量构造预测参数矩阵。步骤4.1中各个imf分量的volterra预测系数的计算方法为：对每个imf分量c重构相空间，令c(n)＝[c(n),c(n-τ),...,c(n-(m-1)t)]为系统的输入向量，c(n+1)为输出向量，则系统的volterra级数展开为：式中，式中vkk(i1,...,ikk)是kk阶volterra核；p是volterra展开阶数，v0为一常系数，则当前volterra展开式为：令volterra级数的预测系数矢量为w(n)，输入信号矢量为z(n)，有：w(n)＝[v0,v1(0),v1(1),...,v1(m-1),v2(0,0),v2(0,1),...,v2(m-1,m-1)]tz(n)＝[1,c(n),c(n-τ),...,c(n-(m-1)τ),c2(n),c(n)c(n-τ),...,c2(n-(m-1)τ)]tc(n+1)＝zt(n)w(n)然后，采用归一化最小均方自适应法求解可得预测系数矢量w(n)的值，其中，m为嵌入维数和τ为延迟时间。本发明的有益效果是：(1)本发明是建立裂纹转子-轴承系统动力学模型，用于研究不同裂纹角下的系统响应，该方法提供了一种严格数学推导的转子-轴承系统的全新的系统响应研究分析方法，可以分析轴承处在非线性油膜力下不同裂纹角对系统响应的影响。(2)在实际应用中难免会受噪声影响，本发明将所得系统响应加入随机白噪声模拟实际信号，采用emd-ciit降噪处理后建立不同裂纹角下的emd-volterra模型级数矩阵，得到矩阵的奇异值作为故障特征量，并加入多组不同白噪声进行数值仿真以验证该方法的稳定性与敏感性，为低转速小裂纹下裂纹故障特征量的提取提供一种思路。附图说明图1是本发明一种转子-轴承系统的横向裂纹故障特征量提取方法的流程图；图2是本发明裂纹转子-轴承系统动力学模型图；图3是本发明裂纹转子示意图；图4是本发明非线性油膜力示意图；图5是本发明产生的系统响应图；图6是本发明给系统响应加入随机白噪声模拟实际信号图；图7是本发明降噪后信号图。具体实施方式本发明一种转子-轴承系统的横向裂纹故障特征量提取方法，流程如图1所示，具体按照以下步骤实施：步骤一：建立裂纹转子-轴承系统动力学模型，获的系统动力学方程，根据系统动力学方程计算不同裂纹角下的系统响应信号；步骤一中所述的裂纹转子-轴承系统动力学模型如图2所示，包括一转轴以及在转轴左右两端呈对称设置的滑动轴承及其转动设置在转轴中间的质量偏心转子，不计圆盘厚度，该裂纹转子-轴承系统所产生的裂纹位于转轴与转子转盘交界处的转轴上。记位于转轴左端滑动轴承的几何中心为o，集中质量为m1，左端滑动轴承在转轴上的中心为o1，x1、y1为左端滑动轴承在转轴x，y方向径向位移，轴承长度l，转子几何中心为o2，集中质量为m2，x2，y2为转子在径向x，y方向位移，e为转子偏心距，滑动轴承的轴套与转轴之间由一层油膜隔开，随着转轴的自转与涡动润滑油在转轴间隙中形成流体动压润滑，油膜厚度为h，r为转轴半径，则裂纹转子-轴承系统动力学方程为：其中，k为转子刚度，c1滑动轴承阻尼系数，c2转子阻尼系数，fx、fy分别为滑动轴承在x、y方向所受到的非线性油膜力，g为重力加速度，β为转子偏心相位，ω为转子转速，t为时间。令x2＝x2/c，y2＝y2/c，x1＝x1/c，y1＝y1/c，τ＝(ωt+β)，ε＝e/c，其中，c为滑动轴承的径向间隙，得系统无量纲化方程为：其中，fx和fyfy分别为滑动轴承在x、y方向所受到的无量纲油膜力。如图3所示，转子刚度在固定坐标系o-x-y下的表达式为：在o'-ζ-η坐标系下，裂纹转子刚度为：式中，k0为完整转子刚度，以裂纹截面圆心o'为原点，垂直横向裂纹的直线为ζ轴，平行横向裂纹的直线为η轴建立旋转坐标系o'-ζ-η，θ为旋转角度，σ为偏心量与旋转坐标轴ζ的夹角，α为裂纹角的一半，δkξ、δkη为ζ、η方向的刚度减少量，式中当转子裂纹深度较浅时，δkξ远大于δkη，可忽略δkη，故为开关函数，通常情况下卧轴转子以在重力作用下的静挠度为主，因此采用重力占优假设符合实际运行情况。用开关函数描述裂纹开闭过程是常用的方法，在保证运算速度和贴近裂纹开闭实际形态的情况下，本发明采用综合模型计算裂纹转子刚度，表达式如下薄油膜位于轴与轴颈之间，如图4所示，由滑动轴承产生的油膜力具有强烈的非线性，可利用经典的capone模型表示油膜力的理论解，将图4中的δ定义为油膜力的起始角位移，δ可表示为：无量纲油膜力fx和fy表示为：式中：步骤二：对不同裂纹角下的系统响应加入白噪声模拟实际信号，并对加入白噪声后的信号进行emd-ciit降噪处理，得到降噪后的信号；在实际环境环境中，噪声的影响难以避免，因此，本文对原始信号加入白噪声模拟实际信号，含白噪声信号图如图6所示，并进行降噪处理，降噪后信号图如图7所示，降噪采用emd-ciit方法；emd方法是希尔伯特-黄(hht)方法的一部分，其关键在于将复杂的非平稳非线性信号分解为一系列内在模态函数，即本征模态函数(imf)；本征模态函数是满足特定条件的函数：(1)在完整信号长度上过零点和极值点的数量相等或不多于一个；(2)任意时刻信号上下包络线关于时间轴对称；一个本征模态函数只有一个波动模式，有研究表明其分解出的imf分量能准确的反映系统的内在特征；步骤二的具体过程为：步骤2.1，对不同裂纹角下的系统响应信号加入白噪声模拟实际信号，得到带噪声信号x(t)；步骤2.2，首先确定x(t)所有极值点，用三次样条插值函数分别连接所有极大值、极小值形成上下包络线，x(t)与上下包络线的均值m1的差记为h1，h1＝x(t)-m1h1；步骤2.3，计算x(t)的第一个imf分量c1(t)：若h1满足imf条件，则h1即为x(t)的第一个imf分量c1(t)；若h1不满足imf条件，则将h1视为原始带噪声信号，重复步骤2.2循环q次，直到h1q＝h1(q-1)-m1q满足imf条件，h1q为第q次时的原始带噪声信号，h1(q-1)为第q-1次时的原始带噪声信号，m1q为h1(q-1)与其上下包络线的均值，此时，可得c1(t)＝h1q，c1(t)为x(t)的第一个imf分量；步骤2.4，从x(t)分离出c1(t)，即就是：r1＝x(t)-c1(t)，将r1当成原始带噪声信号重复步骤2.2-步骤2.3得到x(t)的第二个imf分量c2(t)；步骤2.5，将r2当成原始带噪声信号重复步骤2.2-步骤2.3得到x(t)的第三个imf分量c3(t)，其中，r2＝x(t)-c1(t)-c2(t)，依次类推，直到rl为一单调函数不可再分为止，得到l个本征模态函数由此可知，降噪后信号分解得到一系列本征模态函数ci(t)和余项rl，其中ci(t)包含了频率从高到低的信号成分；步骤2.6，通过小波阈值降噪方法将噪声信号和有用信号分离，获得第一层imf分量中的降噪信号c～1(t)，并计算c1(t)中的真实噪声信号步骤2.7，对噪声信号c1(n)的采样点进行重新排序，得到新的噪声信号c1(a)(t)和新的原始信号其中c1(a)(t)＝alter(c1(n))，alter为选择函数，随机循环和随机置换是最常用的两种选择策略，可用来产生新的噪声序列，它们的区别在于信号样本点的顺序是循环改变还是随机改变。试验结果表明随机循环比随机置换有更好的性能；步骤2.8，对新的原始信号进行emd分解，并对各imf分量采用emd-it降噪处理，得到降噪处理后的各imf分量如下所示：对(表示零点个数)，表示第i个imf中过零瞬间和之间的全部采样点，表示此间隔之间的单一极点。当两个过零点和之间的数据只要有一个极值点超过阈值时，采用重构这之间的数据；如果点和之间的所有极值点小于阈值时这一段全部清零；由此得到emd-it降噪处理后的需要说明的是，这里的ti是一种自适应阈值，表明其可以随着imf的变化而改变；步骤2.9，获得的各个分量进行emd逆过程，即每个分量相加求和即可得到首个降噪信号步骤2.10，然后重复步骤2.7-2.9共z-1次，共获得z个降噪信号，分别为步骤2.11，将步骤2.10求得的所有降噪信号进行求均值，得到最终的降噪后的信号步骤三：对降噪后的信号进行emd分解，提取imf分量；步骤四：对步骤三提取的每个imf分量构造volterra级数模型，得出各imf分量的volterra模型参数后，将各imf分量的volterra模型参数作为矩阵的横矢量构造预测参数矩阵；步骤四的具体过程为：步骤4.1，计算步骤三提取的各个imf分量的volterra预测系数；步骤4.2，得出各imf分量的volterra模型参数后，将其作为矩阵的横矢量构造预测参数矩阵。步骤4.1中各个imf分量的volterra预测系数的计算方法为：对每个imf分量c重构相空间，令c(n)＝[c(n),c(n-τ),...,c(n-(m-1)t)]为系统的输入向量，c(n+1)为输出向量，则系统的volterra级数展开为：式中，式中vkk(i1,...,ikk)是kk阶volterra核；p是volterra展开阶数，v0为一常系数，则当前volterra展开式为：令volterra级数的预测系数矢量为w(n)，输入信号矢量为z(n)，有：w(n)＝[v0,v1(0),v1(1),...,v1(m-1),v2(0,0),v2(0,1),...,v2(m-1,m-1)]tz(n)＝[1,c(n),c(n-τ),...,c(n-(m-1)τ),c2(n),c(n)c(n-τ),...,c2(n-(m-1)τ)]tc(n+1)＝zt(n)w(n)然后，采用归一化最小均方自适应法求解可得预测系数矢量w(n)的值，其中，m为嵌入维数和τ为延迟时间。步骤五：求解所得矩阵的奇异值作为系统的故障特征量。实施例：实施例：在系统的表达式中，通过拟定不同裂纹角，以及设定不同的白噪声，可以验证该方法的有效性。选取：裂纹角α的变化范围取为白噪声信噪比为20db；模型具体参数见表1，表2为1/3倍临界转速裂纹故障特征量，表3为1/2倍临界转速裂纹故障特征量，图5为模型输出系统响应，图6为系统响应加入白噪声信号，图7为降噪后信号。参数数值mb,mr(kg)4,32r,l,c(mm)25,12,0.11cb,cr(ns/m)1050,2100k(n/m)2.5×107e(mm)0.05表1表2表3。当前第1页12