CN112749688A - 一种基于mc罚函数的冲击载荷稀疏识别方法及系统 - Google Patents

Abstract

本发明涉及一种基于MC罚函数的冲击载荷稀疏识别方法及系统，其包括：采集含有噪声信号的转子系统的转子、轴承等振动信号；将振动信号压缩处理后传输至预先建立的稀疏识别模型；对稀疏识别模型进行求解，实现对振动冲击载荷的识别。本发明采用的基于MC罚函数的稀疏分解结果优于L1正则化稀疏分解结果，且在冲击载荷非加载区噪声被很好的抑制。本发明可以广泛在转子系统冲击载荷识别技术领域中应用。

Description

一种基于MC罚函数的冲击载荷稀疏识别方法及系统

技术领域

本发明涉及一种转子系统冲击载荷识别技术领域，特别是关于一种基于MC罚函数的冲击载荷稀疏识别方法及系统。

背景技术

转子系统是大型旋转机械最为重要的组成部分，其可靠性和安全性愈发得到重视，在复杂非平稳的工况下运行时，转子系统的内部结构很容易发生碰撞，引起轴承周期性的振动载荷冲击，导致轴承偏离理想的运动轨迹，并且产生巨大的危害。而其冲击载荷信号很容易被复杂工况下的背景噪声覆盖，再加上采集振动的传感器通常安装在转子系统外壁上，由于复杂的传递路径，信号从产生冲击的轴承到传感器过程中，难免会遇到信号削弱，噪声干扰等状况，因此找到一种合适的燃气轮机转子系统轴承冲击载荷识别方法，对其优化设计及稳定性分析都有很大的帮助。

经过几十年的不断发展，冲击载荷的识别方法已经有了长足的进度，其大致可以分为频域法和时域法，前者发展时间长，主要通过对关心频段上的系统频响函数矩阵求逆来实现冲击载荷识别。比如，Zhou等人基于传递函数的相干分析定义了传递函数相干因子，提出以最小传递函数相干因子为目标的动态载荷识别。时域法恰恰相反，它是在时域内利用阶跃力积分来进行动载荷识别。例如，Chen等人提出了一种时域内基于最优输出跟踪的多源动态冲击载荷识别的方法，通过设计最优输出跟踪器，将载荷识别问题转变为最优输出跟踪问题。Xiao等人提出利用联合奇异熵去噪修正和正则化预优的共轭梯度迭代识别方法，在改善反问题的非适定性的同时也改善冲击载荷识别反问题的求解。Qiao等人提出L1范数的稀疏识别正则化模型，突破基于传统正则化的冲击识别精度不高的瓶颈，Wu等人提出一种基于迭代加权L1正则化稀疏分解方法处理含噪高光谱信息。不仅如此，许多学者尝试在冲击载荷识别方法上添加新理论和新技术，也取得了不少成果，例如最优品质因子、包络解调、小波基、相关峭度解卷积、贝叶斯理论等。然而，上述基于L1范数正则化方法在冲击载荷识别中应用比较广泛，但在高速旋转机械的碰撞冲击中依然存在对响应噪声，高振幅幅值成份上识别精度不高的问题。可见，传统正则化方法识别冲击信号的效率有限，不能满足实时性要求。

发明内容

针对上述问题，本发明的目的是提供一种基于MC罚函数的冲击载荷稀疏识别方法及系统，其稀疏分解结果优于L1正则化稀疏分解结果，且在冲击载荷非加载区噪声被很好的抑制。

为实现上述目的，本发明采取以下技术方案：一种基于MC罚函数的冲击载荷稀疏识别方法，其包括：步骤1、采集含有噪声信号的转子系统的转子、轴承等振动信号；步骤2、将振动信号压缩处理后传输至预先建立的稀疏识别模型；步骤3、对稀疏识别模型进行求解，实现对振动冲击载荷的识别。

进一步，所述步骤2中，采用基于MC罚函数的稀疏识别模型。

进一步，所述基于MC罚函数的稀疏识别模型的建议方法包括以下步骤：

步骤21、确定Tikhonov正则化的目标函数；

步骤22、采用MC罚函数作为目标函数的惩罚项，对目标函数进行正则化；

步骤23、构建基于MC罚函数的稀疏识别模型的目标函数。

进一步，所述MC罚函数为：

式中，a表示MC罚函数的近端算子参数，f表示待识别冲击载荷矢量。

进一步，所述基于MC罚函数的稀疏识别模型的目标函数为：

式中，m表示采集信号的总采样点数，i表示采样点个数，φ(f:a)表示MC罚函数，λ表示正则化参数的加权因子，f表示待识别冲击载荷矢量，y表示稀疏识别结果矢量，H表示函数传递矩阵，a表示MC罚函数的近端算子参数。

进一步，所述步骤3中，采用改进原对偶内点法对稀疏识别模型进行求解。

进一步，采用所述改进原对偶内点法对稀疏识别模型进行求解的方法包括以下步骤：

步骤31、根据拉格朗日的对偶原则，构建作为对偶可行性变量的对偶间隙；

步骤32、引入内点法求解对偶间隙η中约束条件的函数ψ(f,φ)；

步骤33、根据约束条件的函数ψ(f,φ)计算出中心路径轨迹。

进一步，所述步骤32中，函数ψ(f,φ)为：

式中，u表示松弛变量,

进一步，所述步骤33中，计算公式为：

其中，t为势垒参数。

一种基于MC罚函数的冲击载荷稀疏识别系统，其包括：信号采集模块、信号处理模块和模型求解模块；所述信号采集模块，采集含有噪声信号的转子系统的转子、轴承等振动信号；所述信号处理模块，将振动信号压缩处理后传输至预先建立的稀疏识别模型；所述模型求解模块，对稀疏识别模型进行求解，实现对振动冲击载荷的识别。

本发明由于采取以上技术方案，其具有以下优点：本发明利用MC罚函数优化Tikhonov正则化模型，建立适用燃汽轮机转子系统轴承故障的凸优化稀疏分解模型，随后用该模型对仿真信号进行分析，基于MC罚函数的稀疏分解结果优于L1正则化稀疏分解结果，且在冲击载荷非加载区噪声被很好的抑制。

附图说明

图1是本发明具体实施方式中基于MC罚函数的稀疏识别方法流程示意图；

图2是本发明具体实施方式中在a＝1/2时，L1范数与MC罚函数曲线对比图；

图3是本发明具体实施方式中不同取值a的MC罚函数曲线对比图；

图4是本发明具体实施方式中稀疏分解结果时域图；

图5是本发明具体实施方式中稀疏分解过程迭代结果；

图6a是本发明具体实施方式中稀疏分解幅值谱图；

图6b是本发明具体实施方式中稀疏分解包络谱图；

图7是本发明实施例1中试验轴承振动信号时域图；

图8是本发明实施例1中稀疏分解结果时域图；

图9a是本发明实施例1中稀疏分解幅值图；

图9b是本发明实施例1中稀疏分解包络谱图；

图10是本发明实施例2中原始振动信号预处理后时域图；

图11a是本发明实施例2中MC罚函数稀疏分解时域图；

图11b是本发明实施例2中MC罚函数稀疏分解包络谱图。

具体实施方式

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例的附图，对本发明实施例的技术方案进行清楚、完整地描述。显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于所描述的本发明的实施例，本领域普通技术人员所获得的所有其他实施例，都属于本发明保护的范围。

转子系统是大型旋转机械的重要的组成部分，转子、轴承等作为转子系统的关键部件，对转子系统的健康运行起着至关重要的作用，轴承在高温、高压、高转速等恶劣条件下运行，工况环境复杂，极容易发生故障，特别是当轴承受到冲击载荷影响时，严重影响轴承的剩余使用寿命，因此受冲击而造成的轴承摩擦力受力不均是轴承提前失效的主要因素之一。本发明提出一种MC罚函数的冲击载荷稀疏识别方法稀疏分解的冲击载荷识别方法(其中，MC的全称为Minimax Concave penaltyfunction，极大极小罚函数)，该方法利用MC罚函数重新构建正则化稀疏识别模型，然后利用改进原对偶内点法求解该模型，实现振动冲击载荷的识别，实验通过与L1正则化形成对比，结果表明：基于MC罚函数的稀疏分解结果优于L1正则化稀疏分解结果，且在冲击载荷非加载区噪声被很好的抑制。该方法可应用于燃气轮机转子振动信号早期故障诊断中。下面结合附图和实施例对本发明进行详细的描述。

如图1所示，本发明提供一种基于MC罚函数的冲击载荷稀疏识别方法，其包括以下步骤：

步骤1、采集含有噪声信号的转子系统的转子、轴承等振动信号；

步骤2、将振动信号压缩处理后传输至预先建立的稀疏识别模型；

在本实施例中，采用基于MC罚函数的稀疏识别模型；

步骤3、对稀疏识别模型进行求解，实现对振动冲击载荷的识别；在本实施例中，采用改进原对偶内点法对稀疏识别模型进行求解。

上述步骤2中，基于MC罚函数的稀疏识别模型是在Tikhonov正则化基础上利用MC罚函数重新构建正则化稀疏识别模型，使正则化方法能有效避免产生局部最小值，更容易获得次全局最优解。

基于MC罚函数的稀疏识别模型的建议方法包括以下步骤：

步骤21、确定Tikhonov正则化的目标函数；

Tikhonov正则化方法通过添加约束条件使不适定问题变得适定，在残差范数||Hf-y||和解范数||f||₁之间寻找一个称作正则化参数的加权因子λ，进而修正奇异值以提升解的稳定性。Tikhonov正则化的目标函数为：

式中，f表示待识别冲击载荷矢量，y表示稀疏识别结果矢量，H表示函数传递矩阵。

由于振动信号中含有噪声，则式(1)为：

式中，δ表示噪声水平。进一步可以写作无约束的等价形式：

式中，||f||₁表示L1范数惩罚项，其在残差

和稀疏正则化算子||f||₁间建立了折中关系。式(3)中的L1范数项||f||₁具有促使待识别冲击载荷中较小数值变成零的作用，从而诱导得到冲击载荷的稀疏解。

步骤22、采用MC罚函数作为目标函数的惩罚项，对目标函数进行正则化；

MC罚函数公式如下：

式中，a表示MC罚函数的近端算子参数。

由于在实际应用中，求解维度较大时，或工况复杂多变时，除了固定的阈值函数外，L1阈值函数与同一函数有一定的距离，这使得L1惩罚函数导致信号幅度被低估，该状况下有一些微弱冲击可能会被忽略，对早期的故障检测容易造成误判。故本实施例中采用MC罚函数，MC罚函数的近端算子与同一函数表示一致，这样就能保证信号幅度不会被低估。

不同近端算子参数的MC罚函数曲线图如图2，图3所示。由图2可以看出，MC罚函数与L1范数的不同之处在于MC罚函数不以x的正负数为判断条件，而是通过引入数值a来作为分割点，且在|x|的基础上添加约束项，其他情况下也不再为零，这样相对于L1范数而言，在迭代过程中不会盲目的将数值化为零，在一定程度上保留了更多信息。从图3中可以得出解论，当a越靠近0时，MC罚函数越接近L1范数，当a接近1时，MC曲线的稳定约束区域越来越大，不难看出，当a＝1或者0.5时，MC罚函数的约束效果最好，但考虑到稳定性以及防止过拟合问题，最终选取a＝0.5作为MC罚函数的近端算子参数。

步骤23、构建基于MC罚函数的稀疏识别模型的目标函数为：

式中，m表示采集信号的总采样点数，i表示采样点个数，φ(f:a)表示MC罚函数，λ表示正则化参数的加权因子，f表示待识别冲击载荷矢量，y表示稀疏识别结果矢量，H表示函数传递矩阵，a表示MC罚函数的近端算子参数。

上述步骤3中，采用改进原对偶内点法对稀疏识别模型进行求解的方法包括以下步骤：

步骤31、根据拉格朗日的对偶原则，构建作为对偶可行性变量的对偶间隙；

MC罚函数均属于非凹罚函数，该类正则化模型不存在显式解而不得不依赖优化算法，故不能用奇异值分解方法求解模型。传统的原对偶内点法在求解模型时存在对噪声敏感以及高振幅幅值被低估现象，所以本发明采用改进原对偶内点法求解冲击载荷的稀疏识别模型，鲁棒性强，收敛迅速。

其中，拉格朗日对偶公式为：

式中，G(v)表示拉格朗日对偶函数，v表示Lagrange乘子，T表示转置矩阵，i＝1,…m。

该问题极小值存在的必要条件是拉格朗日函数对所有变量及乘子的偏导数为0，由于此方程组为线性方程，可用共轭梯度法进行求解。

将基于MC罚函数的稀疏识别模型的目标函数和相应的拉格朗日对偶函数之间的差异称为对偶间隙η，其公式如下：

步骤32、引入内点法求解对偶间隙η中约束条件的函数ψ(f,φ)为：

式中，u表示松弛变量,

步骤33、根据约束条件的函数ψ(f,φ)计算出中心路径轨迹；

计算公式为：

其中，t表示一个小的正标量，称为势垒参数。

在一个优选的实施例中，本发明基于MC罚函数的冲击载荷稀疏识别方法还包括：

步骤4、基于仿真信号对载荷识别方法进行验证；

在本实施例中，采用模拟仿真冲击载荷的振动信号进行验证，为了提高仿真信息的有效性，冲击载荷的仿真信号模型如下：

式中，β表示阻尼系数，ω表示振动频率。其他参数设置采样频率f_s为2000Hz，延迟时间t₀为0.05s,信噪比SNR＝5.7724，经计算，故障特征频率为300Hz。仿真信号的时域图如图2所示。

采用MC罚函数稀疏识别模型对加噪仿真信号进行稀疏分解识别冲击载荷信号，其稀疏分解结果如图4、图5所示。

从图4可以看出，基于MC稀疏分解结果受噪声影响较小，而基于L1正则化分解结果受噪声影响较大，并且存在高幅值在一定程度上被削弱的现象，这主要因为L1范数误差的线性增长速度过快引起的。从收敛速度角度分析，图5中可以看出，虽然基于L1范数正则化和基于MC稀疏分解的迭代过程均趋于收敛，但收敛速度存在明显差距，基于L1范数正则化在迭代30次左右开始趋于稳定，且迭代过程中存在阶梯式收敛，而基于MC范数的稀疏分解在迭代20次左右便趋于稳定，并且在迭代过程中稳定性有明显提升。通过比较MC范数稀疏分解和L1范数正则化的稀疏结果，可以得出基于MC范数稀疏分解方法不仅能降低噪声对信号的干扰，而且还能保留了部分高振幅成分。

图6a是基于L1范数和MC罚函数正则化稀疏分解的幅值图，通过对比两种方法，能清楚看到故障特征频率为300Hz，与模拟信号故障特征频率300Hz相等，除此之外，从图中也能看到噪声因素对识别后结果的影响程度，基于MC罚函数的技术分解较好的抵挡噪声对信号的干扰，且高振幅幅值保留较好，而基于L1范数正则化方法在噪声因素的敏感度上较强，且在高振幅上有明显的削弱，为了能更清楚的看出两者的差距，继续进行包络谱分析，图6b是基于L1范数和MC罚函数的正则化包络谱图，从结果上来看，基于MC罚函数的稀疏分解在噪声敏感度上较弱，高振幅上有一定保留。

本发明还提供一种基于MC罚函数的冲击载荷稀疏识别系统，其包括信号采集模块、信号处理模块和模型求解模块；

信号采集模块，采集含有噪声信号的转子系统的转子、轴承等振动信号；

信号处理模块，将振动信号压缩处理后传输至预先建立的稀疏识别模型；

模型求解模块，对稀疏识别模型进行求解，实现对振动冲击载荷的识别。

实施例1：

采用单次振动冲击实验验证本发明识别方法的有效性和适应性。

本次实验数据集采集来自于电主轴试验台的数据。电主轴试验台由电主轴、液压加载装置、测功机、油气润滑装置组成。采集放置于主轴模拟刀柄位置(加载棒)加速度传感器信息，获取主轴径向和轴向位置信号。由于试验台的抱夹机构，使得径向跳动小于轴向跳动，所以分析时采用轴向方向的振动信号。

对主轴施加不同瞬间载荷来模拟主轴轴承的振动冲击情况，试验中采样频率f_s为2000Hz，采集工况为1000空载，主轴转速恒为2000r/min，振动信号采用东华测试采集分析仪记录，采集通道选用1通道，试验台数据结果的时域信号如图7所示。

从图7时域图中可以看出时域波形杂乱且含有大量的噪声成分，振动冲击脉冲不明显，为了有效提取轴承振动冲击信号特征及验证该方法的有效性，分别采用L1正则化稀疏分解法和MC罚函数稀疏分解算法进行稀疏分解，结果如图8所示。

从图8两种方法处理后时域图中可以看出，MC罚函数稀疏分解在高振幅成份上保留的效果优于L1正则化，分别为1.8μm和0.96μm。噪声控制范围也小于L1正则化。接下来分别对以上两种结果进行幅值和包络谱分析，结果如图9a、图9b所示。

从图9a中可以清楚的看到振动冲击载荷信号均被保留下来，故障特征频率都出现在320Hz附近，但同时可以发现两者在振幅和噪声方面存在较大差异，基于MC罚函数算法高振幅保留到1.1，而基于L1正则化算法保振幅仅保留到0.27，两者相差明显，噪声方面，L1正则化幅值图中在400Hz和900Hz附近出现幅值波动，不符合本实验单次振动冲击载荷的试验预期，而MC罚函数的幅值图中却没有出现这中现象，噪声控制效果明显提升。

包络谱分析的图9b可反映出两者的差异：L1正则化幅值图幅值为0.007，噪声方面在670Hz附近出现明显振幅波动；MC罚函数幅值图幅值为4.1，其他频率无明显振幅波动。

通过本次实验验证，MC罚函数稀疏分解算法能很好的控制噪声对振动冲击信号的干扰以及在高振幅幅值保留成分上的优势，证明了该方法的有效性和适应性。

实施例2：

采用多次振动冲击实验验证本发明识别方法的有效性和适应性。

为了进一步验证该方法的适应性，用某型号转子系统试验台的数据，其基本参数为：工作转速6000r/min，传感器测试记录高压轴的转频，采样频率f_s为6000Hz，运行速度保持在0.1工况下运行，测试过程中，从低转速缓慢升至0.1工况，其传感器检测到后侧点转子轴承出现振幅明显增大的现象。振动信号利用亿恒的INTEDER数据采集与分析仪记录，采集通道选用16通道，采集设备采集到的信号经数据压缩以及剖析处理后的时域图如图10所示。

基于L1惩罚函数和MC惩罚函数提取的振动冲击时域和包络谱分析如图11a、图11b所示。从图11a中对比可以看出，在多次振动冲击状况下，噪声对两种方法的适应性差距非常明显，振动冲击现象主要发生在0.1s、0.24s、0.49s、0.87s四个时间点附近，而从结果来看，基于L1正则化稀疏分解结果除以上四个时间点外均出现明显的异常噪声干扰性幅值，冲击识别不明显，处理效果不理想。而基于MC罚函数的稀疏分解结果没有明显的异常噪声干扰，噪声得到很好的控制，0.1s、0.24s、0.49s、0.87s四个时间点的振动冲击识别效果明显，且高振幅的幅值保留成分明显优于L1正则化方法。不仅如此，如图11b所示，从包络谱分析中两者的差异也能很好的看出。

通过本次实验验证，在多次振动冲击的情况下，MC罚函数稀疏分解算法控制噪声对冲击信号的干扰能力优势更加突出，噪声控制在合理的范围内，对稀疏分解过程的干扰性低，其二，在高振幅幅值保留成分上优势更加明显，该实验进一步证明了该方法的有效性和适应性。

综上所述，本发明不仅可以分离出轴承故障信号中的冲击成分，而且有效抑制噪声对稀疏分解过程的干扰，具有良好的降噪效果，同时也更好的保留了高振幅幅值成分。而且本发明采用的改进后的原对偶内点法提高了稀疏识别模型迭代收敛速度，且过程收敛平缓，与L1正则化方法对比可知，本发明鲁棒性好，训练过程更加稳定。

本领域内的技术人员应明白，本申请的实施例可提供为方法、系统、或计算机程序产品。因此，本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且，本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。

本申请是参照根据本申请实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的装置。

这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。

这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。

Claims (10)

1.一种基于MC罚函数的冲击载荷稀疏识别方法，其特征在于，包括：

步骤1、采集含有噪声信号的转子系统的转子、轴承等振动信号；

步骤2、将振动信号压缩处理后传输至预先建立的稀疏识别模型；

步骤3、对稀疏识别模型进行求解，实现对振动冲击载荷的识别。

2.如权利要求1所述识别方法，其特征在于，所述步骤2中，采用基于MC罚函数的稀疏识别模型。

3.如权利要求2所述识别方法，其特征在于，所述基于MC罚函数的稀疏识别模型的建议方法包括以下步骤：

步骤21、确定Tikhonov正则化的目标函数；

步骤22、采用MC罚函数作为目标函数的惩罚项，对目标函数进行正则化；

步骤23、构建基于MC罚函数的稀疏识别模型的目标函数。

4.如权利要求3所述识别方法，其特征在于，所述MC罚函数为：

式中，a表示MC罚函数的近端算子参数，f表示待识别冲击载荷矢量。

5.如权利要求3所述识别方法，其特征在于，所述基于MC罚函数的稀疏识别模型的目标函数为：

式中，m表示采集信号的总采样点数，i表示采样点个数，φ(f：a)表示MC罚函数，λ表示正则化参数的加权因子，f表示待识别冲击载荷矢量，y表示稀疏识别结果矢量，H表示函数传递矩阵，a表示MC罚函数的近端算子参数。

6.如权利要求1至5任一项所述识别方法，其特征在于，所述步骤3中，采用改进原对偶内点法对稀疏识别模型进行求解。

7.如权利要求6所述识别方法，其特征在于，采用所述改进原对偶内点法对稀疏识别模型进行求解的方法包括以下步骤：

步骤31、根据拉格朗日的对偶原则，构建作为对偶可行性变量的对偶间隙；

步骤32、引入内点法求解对偶间隙η中约束条件的函数ψ(f，φ)；

步骤33、根据约束条件的函数ψ(f，φ)计算出中心路径轨迹。

8.如权利要求7所述识别方法，其特征在于，所述步骤32中，函数ψ(f，φ)为：

式中，u表示松弛变量，

9.如权利要求7所述识别方法，其特征在于，所述步骤33中，计算公式为：

其中，t为势垒参数。

10.一种基于MC罚函数的冲击载荷稀疏识别系统，其特征在于，包括：信号采集模块、信号处理模块和模型求解模块；

所述信号采集模块，采集含有噪声信号的转子系统的转子、轴承等振动信号；

所述信号处理模块，将振动信号压缩处理后传输至预先建立的稀疏识别模型；

所述模型求解模块，对稀疏识别模型进行求解，实现对振动冲击载荷的识别。

Images