基于IMM-KEnPF的风机桨距系统故障诊断方法

技术领域

本发明涉及故障诊断技术领域，特别是涉及基于IMM-KEnPF的风机桨距系统故障诊断方法。

背景技术

风电机组中的变桨距系统是故障常发部件之一，故障会影响风机功率输出的稳定性，并会导致转子、转轴和轴承上的负载不平衡，给风机的安全运行带来隐患。因此，对风机变桨距系统的故障诊断方法进行研究，对于保障风电机组的安全稳定运行具有重要作用。

相对于风电机组的其他子系统，变桨距系统故障诊断方法的研究成果较少。ZhangYoumin最先提出了基于IMM的故障诊断的基本框架，该方法对于系统参数或者结构会发生变化的系统也可以较好地完成自适应估计。胡艳艳等提出了网络化异步IMM融合滤波算法，将未知的故障幅值作为增广的系统状态变量，根据模型概率检测和定位故障，同时得到故障幅值和系统状态的联合估计。王剑等将无迹卡尔曼滤波与交互式多模型相结合，提出了IMM-UKF故障诊断算法。GADSDEN S A等使用时变平滑边界层(VBL)来改进平滑变结构滤波器(SVSF)，并与交互式多模型结合设计出基于IMM-SVVF-VBL的故障诊断算法。但是UKF和SVSF处理非线性的能力有限，诊断性能有待进一步提升。Zhengjiang提出将PF与交互式多模型方法结合，设计出一种适用于强非线性系统的故障诊断方法。但PF中存在的粒子退化和样本匮乏导致的估计精度下降以及计算耗时问题。同时，交互式多模型需要根据所建立的模型集，集成多个滤波器，模型集中的子模型越多所需要的滤波器就越多，计算负荷会急剧增加。

针对以上问题，目前的研究中主要使用群智能算法对粒子后验分布进行优化，但是群智能优化容易陷入局部最优，且寻优迭代过程会耗费许多计算时间，而使用次优滤波器直接优化建议分布改善重要性采样过程则更直观和高效。其中常见的次优滤波器是无迹卡尔曼滤波(UKF)。Merwe等将无迹卡尔曼滤波(Unscented Kalman Filter,UKF)引入粒子滤波，根据最新的观测值通过无迹卡尔曼滤波来预测状态值的变化，以此来优化粒子滤波的建议分布，缓解粒子退化现象。但是UKF中的Sigma点的计算时间会随着维度和粒子增加而急剧增加。Delft G V等设计出一种选取随机点，用来预测状态后验分布作为粒子滤波的建议分布，称作集合粒子滤波(Ensemble Particle Filter,EnPF)，有着更好的处理系统非线性的能力和计算效率。当粒子数量增加时，使用次优滤波器优化粒子滤波的方法同时也增加了算法的计算复杂度，会严重影响算法的实时性，如何在提高粒子滤波算法估计精度的同时从算法层面提高算法的实时性仍是一个研究难题。

发明内容

本发明提供的基于IMM-KEnPF的风机桨距系统故障诊断方法，可以解决现有技术中存在的问题。

本发明提供了基于IMM-KEnPF的风机桨距系统故障诊断方法，包括以下步骤：

步骤1，输入交互；风机变桨距系统的多种故障模型构成故障模型集，故障模型集中的故障模型之间随机跳变切换，前一时刻的故障模型到下一时刻的故障模型之间的切换服从故障模型转移概率；基于交互式多模型IMM，利用上一时刻的各故障模型概率和故障模型转移概率计算故障模型之间的交互概率，再使用交互概率融合各粒子滤波器的上一时刻的滤波信息，得到各故障模型在上一时刻初始的估计交互值和对应的误差协方差矩阵；

步骤2，并行滤波；利用步骤1得到的各故障模型在上一时刻初始的估计交互值和对应的误差协方差矩阵，为各个粒子滤波器生成粒子集，按照集合卡尔曼滤波EnKF优化的建议分布以及变桨距系统的观测值，对粒子集中的粒子进行重要性采样，得到粒子的预测集合；在重要性采样过程中，利用库尔贝克-莱布勒散度KLD采样算法自适应调整采样粒子的数量；

得到预测集合后，根据变桨距系统的观测值计算粒子权值，对粒子残值利用残差重采样方法进行重采样，得到新的粒子集合，将新的粒子集合中粒子的粒子均值、各个粒子滤波器状态估计的误差协方差矩阵、估计残差和残差协方差矩阵作为输出；

步骤3，模型概率更新；根据步骤2输出的估计残差和残差协方差矩阵建立各故障模型滤波输出的似然函数，使用该似然函数计算各个故障模型的模型概率；

步骤4，融合估计；得到各个故障模型的模型概率后，将其与步骤2输出的粒子均值和各个粒子滤波器状态估计的误差协方差矩阵进行融合，得到估计结果；

步骤5，故障诊断；基于步骤3得到的各个故障模型的模型概率，利用决策函数实现变桨距系统的故障诊断。

本发明中的基于IMM-KEnPF的风机桨距系统故障诊断方法，针对变桨距系统的非线性非高斯的特点引入粒子滤波算法，然后从准确性和实时性两个方面对其进行了改进，最后将交互式多模型方法和粒子滤波算法结合。首先引入集合卡尔曼滤波(EnKF)优化PF的建议分布，缓解PF的粒子退化问题，提高PF的估计精度；其次引入KLD采样，根据后验估计分布与真实分布的KLD来自适应的调整所需要采样的粒子数目，提高PF的计算效率和实时性。最后将改进的PF与IMM相结合，在风电机组仿真平台上实现了变桨距系统多传感器故障的诊断。

附图说明

为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明中诊断方法的流程图；

图2为仿真实验框架；

图3为模型概率的变化曲线；

图4为故障诊断结果示意图；

图5为状态估计结果示意图；

图6为本发明方法和对比方法的计算时间对比。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

在介绍本发明的方法前，首先对风机变桨距系统的故障进行介绍。

将风机变桨距系统考虑成一个活塞伺服系统，进而将其建模成一个二阶传递函数模型：

式中，ζ为变桨距系统的阻尼系数，ω

根据上述的基准模型对不同的故障建立相应的数学模型。将变桨距系统的数学模型从传递函数的形式转化成状态空间的形式：

式中，f(·)为变桨距系统的状态转移函数，h(·)为观测函数，w

表1传感器故障模型

参照图1，本发明提供了基于IMM-KEnPF的风机桨距系统故障诊断方法，该方法包括以下五个步骤：

步骤1，输入交互。风机变桨距系统的多种故障模型构成故障模型集，故障模型集中的故障模型之间是随机跳变切换的，前一时刻的故障模型到下一时刻的故障模型之间的切换服从故障模型转移概率。基于交互式多模型IMM，利用上一时刻的各故障模型概率和故障模型转移概率计算故障模型之间的交互概率，再使用交互概率融合各粒子滤波器的上一时刻的滤波信息，得到各故障模型在上一时刻初始的估计交互值和对应的误差协方差矩阵，实现粒子滤波器的初始化。

步骤2，并行滤波。利用步骤1得到的各故障模型在上一时刻初始的估计交互值和对应的误差协方差矩阵，为各个粒子滤波器生成粒子集，按照集合卡尔曼滤波EnKF优化的建议分布以及变桨距系统的观测值，对粒子集中的粒子进行重要性采样，得到粒子的预测集合。在重要性采样过程中，利用库尔贝克-莱布勒散度KLD采样算法自适应调整采样粒子的数量。

得到预测集合后，根据变桨距系统的观测值计算粒子权值，对粒子残值利用残差重采样方法进行重采样，得到新的粒子集合，将新的粒子集合中粒子的均值、各个粒子滤波器状态估计的误差协方差矩阵、估计残差和残差协方差矩阵作为输出。

步骤3，模型概率更新。根据步骤2输出的估计残差和残差协方差矩阵建立各故障模型滤波输出的似然函数，使用该似然函数计算各个故障模型的模型概率。

步骤4，融合估计。得到各个故障模型的模型概率后，将其与步骤2输出的粒子均值和各个粒子滤波器状态估计的误差协方差矩阵进行融合，得到融合估计结果。

步骤5，故障诊断。基于步骤3得到的各个故障模型的模型概率，利用决策函数实现变桨距系统的工作模式诊断。

下面对每个步骤进行详细说明。

步骤1，由多个故障模型建立的故障模型集为M＝{M

式中：M

P{M

式中，π

将交互式多模型IMM与粒子滤波KEnPF进行结合时，交互环节的实施有模型层面的、粒子层面的以及前两者相结合的三种方式。现有文献对以上三种交互的实施方式做了对比实验，结果表明第一种方式的目标跟踪的性能最好，因此，本发明在故障诊断设计上也采取该方式。

在步骤1迭代周期开始时，利用上一时刻的各故障模型概率和故障模型转移概率计算故障模型之间的交互概率，再使用交互概率融合各粒子滤波器的上一时刻的滤波信息进行融合，实现滤波器的初始化，交互的过程如式(5)-(9)所示：

式中：i,j∈M，

步骤2，首先，利用前一步的交互结果

其中，rand表示随机函数。

随后，进入改进粒子滤波KEnPF的流程，详细流程如下：

1、初始化。在k＝0时刻，从先验概率分布中随机采样N个权值大小相同的粒子

2、定义背景集合。为第i个粒子随机采样n个集合点，构成背景集合

3、更新背景集合。

随机采样n个初始样本集合

根据式(12)计算背景集合的均值，并用式(13)和式(14)代替方差计算。

通过下式计算卡尔曼增益：

式中，R

最后根据下式计算分析集合的均值和方差：

使用式(11)-(16)更新背景集合得到分析集合

4、重要性采样和权值更新。从EnKF优化的建议分布，即上述得到的各个粒子滤波器生成粒子集中采集粒子，得到粒子的预测集合：

将标准粒子滤波器中的状态转移函数作为建议分布替换，使得每个粒子

由于充分利用最新时刻的观测信息，获得的建议分布将更合理，可改善粒子滤波的粒子退化现象和样本匮乏，更易于应对系统状态的突变。

对预测集合中的第i个粒子进行重要性采样

S＝S+ω

5、更新自适应粒子数目。判断x

6、当i＜n

7、权值归一化和重采样。使用权值累加和进行权值归一化：

得到粒子集

8、加权输出估计结果：

KLD采样的方法为：

任意两个分布p与q之间的KLD由式(25)定义：

利用KLD调整粒子采样规模，就是要使粒子的采样数目满足真实后验分布与估计后验分布之间的KLD小于ε的概率为1-δ。假设有n个粒子从k个离散的空间中随机地取出，用

其中，

式中，z

根据观测值计算残差，用残差计算粒子权值，并归一化。

式中：z

最后，分别计算各个粒子滤波器状态估计的误差协方差矩阵

步骤3，使用各故障模型滤波输出的似然函数来计算故障模型概率，如式(35)所示。

式中，似然函数根据各滤波器估计结果的残差

故障模型概率在交互式多模型故障诊断算法中起到关键性的作用，除了利用其对滤波器的滤波结果进行融合输出，也需要被用来进行故障模型的匹配进行诊断。在系统处于工作正常的情况下，与正常模式对应故障模型的故障模型概率应该接近于1，当工作模式发生变化，或发生故障时，对应的匹配故障模型的故障模型概率也就进阶地上升并接近于1，而其他故障模型的故障模型概率下降到接近于0。

步骤4，在获得故障模型概率之后，使用故障模型概率对各滤波器的滤波结果进行加权融合，得到IMM-KEnPF的估计结果，如(37)-(39)所示。

步骤5，对于故障模型概率所蕴含的模式信息，可以利用式(40)-(41)的决策函数实现工作模式的识别和诊断。当某一个故障模型与当前系统匹配，而其他故障模型都失配时，说明系统正处于对应的故障模式当中。

式中，μ′

仿真实验

为了运用所提出的基于IMM-KEnPF方法解决变桨距系统的故障诊断问题并验证所提出改进的效果，本发明建立了如图2所示的实验框架，先通过Matlab/Simulation搭建的风机仿真平台采集变桨距系统含有故障的运行数据，再将数据传给交互式多模型进行故障诊断和状态估计。

首先，根据建立的故障模型集，对应设计4个KEnPF滤波器，通过各滤波的滤波结果进行故障模型概率的计算，从而进行变桨系统的故障诊断和状态估计，最后将此方法与基于MMAE-PF以及IMM-PF的故障诊断方法进行了对比。实验参数设置如下：仿真时间为T＝6s，采样间隔ΔT＝0.01，粒子滤波器粒子数N＝50，过程噪声和观测噪声分别为w～Γ(0.1,0.1)，v～N(0,0.001)，KEnPF中最小采样粒子数n

实验硬件为Intel(R)Corei5-8250U CPU@1.6GHz，内存8G，操作系统为Windows10，编程软件为Matlab R2017b。

方法改进前后的故障模型概率变化对比如图3所示。可以看到在200ΔT处故障模型切换后的故障模型概率发生了反转，故障模型0的故障模型概率上升到接近1，其他故障模型的概率下降到接近0。而改进后的故障诊断算法相较标准IMM-PF算法，由于KEnPF改进的估计性能，故障模型概率也更接近于1。由于诊断和融合估计都要使用故障模型概率，这样的效果有利于提高算法诊断的准确度和状态估计的精度。从图4改进后的诊断结果可以看出，正确识别故障的数量明显增加，故障模型失配得到了减少；从图5改进后的状态估计结果可以看出桨距角也能更好地逼近真实值。

图4中，从上到下分别是故障的预设模式序列、基于MMAE-PF的故障诊断结果、基于标准的IMM-PF故障诊断结果、基于本发明的IMM-KEnPF的故障诊断结果。可以看出基于MMAE-PF方法的故障诊断结果相较于其他两种方法有较多的误报和漏报，基于IMM-PF的故障诊断结果中，除了在故障模式发生切换时存在一定的时间延迟，整个过程有较少的漏报，诊断结果整体上与预设的模式序列相吻合。

图5是MMAE-PF、IMM-PF和IMM-KEnPF三种算法的自适应估计结果，从中可以看出MMAE-PF的估计偏离实际桨距角的程度较大，而IMM-PF大部分区间内的桨距角估计较为准确，误差较小，只是在500～600ΔT发生传感器卡死时，桨距角估计都发生了振荡。这是由于输入值在变化，而滤波器无法获得相应的输出值的更新。图6显示了IMM-PF与IMM-KEnPF的整个故障诊断方法计算的时间对比，可以看出相较于固定粒子数的IMM-PF方法，IMM-KEnPF方法由于具有采用KLD采样对粒子数进行削减，而在计算时间上有大幅的节省，并且随着粒子数的提升，效果越明显。

从图3～6的实验结果可以看出，本发明提出的基于IMM-KEnPF的变桨距传感器的多故障诊断方法总体上是有效的，而在此基础上，对交互式多模型和滤波器的改进进一步提升了故障诊断和状态估计性能。为综合评价这两种算法的性能，选取了以下三个指标：

(1)正确诊断(CDID)：诊断结果与实际模式匹配的个数；

(2)均方根误差(RMSE)：如式(42)所示的交互式多模型的融合估计的均方根误差；

式中：x

(3)计算时间(runtime):每个故障诊断算法运行周期的计算时间。本发明选取粒子数为50的情况下统计计算时间。

将改进前和改进后的故障诊断方法各独立运行50次蒙特卡洛实验，计算各次的平均值，得到的各性能指标统计结果如表2所示。

表2 IMM-KEnPF与IMM-PF性能对比

从表中的数据可以看出，改进后的故障诊断方法提高了18.22个正确诊断数，提升了近3％，降低了56.17％的估计误差，并且计算时间减少了72.2％。此统计结果表明本发明所提出的IMM-KEnPF的故障诊断方法从诊断准确性、状态估计精度以及实时性都得到了较好的提升。

尽管已描述了本发明的优选实施例，但本领域内的技术人员一旦得知了基本创造性概念，则可对这些实施例作出另外的变更和修改。所以，所附权利要求意欲解释为包括优选实施例以及落入本发明范围的所有变更和修改。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围。这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。