卢俊杰,黄金泉,鲁 峰
(南京航空航天大学能源与动力学院,南京210016)
涡扇发动机结构复杂,工作环境恶劣,属于1种故障多发系统。据统计,涡扇发动机气路部件故障占涡扇发动机总体故障的90%以上,因此,实时检测发动机健康状况并进行气路性能分析是提高发动机安全可靠性的重要途径[1-3]。发动机气路部件的效率变化系数、流量变化系数等健康参数是发动机气路故障的状态特征,会直接导致转速、温度、压力等发动机测量参数的变化,因此发动机气路部件故障诊断主要采用特定的算法根据测量参数的变化来估计出健康参数[4],从而分析发动机气路部件的健康状况。发动机故障诊断方法主要有基于模型(如最小二乘方法[5],卡尔曼滤波方法[6]),基于数据(如神经网络[7],支持向量机[8])和基于知识(如专家系统)的方法。随着发动机部件级模型的精度以及计算机速度的提高,基于模型的故障诊断方法在工程实际中得到了广泛的研究和应用。
卡尔曼滤波算法受限于线性高斯系统,一些改进的卡尔曼算法[9-10]虽然适用于非线性问题,但依然依赖于高斯噪声的假设[11]。而粒子滤波算法适用于含非高斯噪声的非线性系统,于是有学者提出使用粒子滤波算法进行气路部件健康参数的估计,并取得较好效果[12]。标准粒子滤波用于发动机故障诊断主要存在以下2方面问题:(1)标准粒子滤波直接采用先验分布作为重要性密度函数[13],没有结合最新观测值,导致在突变故障下,需要诊断步数多;(2)发动机健康参数的维度高,滤波器的抽样率低,同时出于实时性考虑,粒子数少,所以有效粒子数少,导致诊断结果噪声水平高,估计精度低。
本文针对粒子滤波应用于突变故障诊断时诊断步数多和诊断结果噪声水平高的问题,通过引入伪协方差与自适应调整似然分布尾部的平坦程度,提出使用伪协方差和自适应似然分布结合的1种改进粒子滤波方法。
粒子滤波是1种基于递推贝叶斯估计和蒙特卡罗积分的统计滤波方法,其主要思想是:首先,由系统状态向量的经验条件分布在状态空间中产生1组随机的粒子;其次,根据新的观测值调整粒子的位置和权值;最后,通过调整后的粒子信息修正最初的条件分布[14]。实质即为根据粒子及其权值组成的离散随机测度来近似状态的后验概率分布,且这种近似会随着粒子的数目趋于无穷多时收敛于真实后验概率分布。
考虑如下的非线性系统
式中:xk、uk、zk分别为k时刻的状态量、控制量和观测量;ωk、vk分别为k时刻的过程噪声和观测噪声;f(·)、h(·)分别为状态转移函数和观测函数。
k时刻目标状态xk的后验概率分布p(xk|z1,…zk)可以用{来描述,其中{是权值分别为{的粒子集,N为采样粒子数。设归一化权值为则k时刻系统状态的后验概率分布近似离散加权为
标准粒子滤波算法的重要性密度函数选择为最容易实现的先验概率密度
由归一化权值加权求和得状态的最优估计值
高斯粒子滤波是标准粒子滤波的1种实现方法,通过高斯密度函数来逼近状态后验概率分布[15],基本思想是利用反映粒子平均数值水平的数学期望和反映粒子离散程度的协方差进行递推来简化粒子集的更新过程。本节在高斯粒子滤波实现方法的基础上,提出1种改进的粒子滤波算法。该方法结合发动机健康参数估计的特点,引入伪协方差的代替协方差,以期望减少突变故障诊断步数;另外通过对似然分布自适应调整,以期望降低诊断结果噪声水平。
在发动机故障诊断中,在无突变故障发生的情况下,粒子加权均值基本稳定,可以近似认为。粒子集的协方差∑k为
在似然分布位于转移先验分布尾部或者观测模型具有很高精度时[16],很多样本由于归一化权重很小而成为无效样本,过低的抽样率会导致粒子滤波失败,增大粒子数是其解决方法之一,但随之而来的巨大计算量使该方法在工程应用中无法实现。
发动机故障诊断需要估计多个健康参数,高维度必然导致似然分布处于转移先验尾部,如图1所示。从图中可见,很少的粒子位于B处,绝大部分粒子位于A处成为无效粒子。同时受实时性限制,粒子总数少,所以滤波器有效粒子数低,诊断结果的噪声水平高。针对该问题,本节通过对似然分布进行自适应调整以提高抽样率,使之尾部更为平坦,增加先验和似然的重叠区,降低诊断结果的噪声水平。改变似然函数的分布就是改变粒子权值的分布,所以只需对标准粒子滤波算法的权值更新作自适应改进。
图1 状态转移函数和似然函数的关系
式中:β为比例因子,根据权值离散水平进行调整;κ为挤压因子。根据权值数值水平进行调整,取β∈(0.6,1.5),κ∈(0.3,0.9)。
(3)计算各粒子的权值
(5)将各粒子的权值归一化
(6)求状态估计值xˆk和粒子集伪协方差
为了表示故障发生后部件性能的变化程度,引入旋转部件的效率变化系数ΔSEi和流量变化系数ΔSWi[17]
式中:ηi、Wi为部件的实际效率和流量;为部件效率和流量的理想值。
健康参数选为风扇、压气机、高压涡轮和低压涡轮的效率和流量的变化系数[18],共8个,定义为
将其增广到发动机状态中,则
式中:NL、NH分别为低、高压转子转速;ΔSE1、ΔSE2、ΔSE3、ΔSE4分别为风扇、压气机、高压涡轮、低压涡轮效率变化系数;ΔSW1、ΔSW2、ΔSW3、ΔSW4分别为风扇、压气机、高压涡轮、低压涡轮流量变化系数;上标T表示将行向量转换为列向量的变换,涉及的控制量有:Wfb为燃油量,A8为尾喷管面积。
基于非线性模型的气路部件故障诊断原理如图2所示。通过发动机的输出值与模型预测值之间的残差结合粒子滤波算法对气路部件健康参数进行估计。发动机传感器参数选取见表1,相应传感器的噪声水平根据文献[12,17]确定。
图2 基于粒子滤波的故障诊断原理
表1 传感器选取及其测量噪声
高斯分布为统计学中常用的理想噪声模型,然而现实应用中噪声通常都是比较复杂的,为了模拟发动机运行时的复杂噪声,本文选取由伽马噪声以及高斯噪声简单混合组成的双模噪声,对发动机健康参数突变故障进行估计,以验证粒子滤波在非高斯噪声下状态参数估计的准确性。为了比较的公平性,标准粒子滤波算法以及改进粒子滤波算法的粒子数均设为50。另外,设~N(μm,)~Γ(μm,),观测噪声,其中 m∈{1,2,…,8}表示对应传感器,和相互独立。根据表1中传感器测量噪声标准差取R值,令r=1,μ=-,σ=,=mmmm,则 E()=0,D()=Rm。过程噪声参考文献[12]选取为对角阵 Q=0.00152×I10×10。
在地面标准状况的最大工作状态下,分别用标准粒子滤波器和改进粒子滤波器对发动机气路部件的单部件故障模式和多部件故障模式进行仿真估计。参考NASA在MAPSS仿真平台中涡扇发动机完成一定工作循环数量后气路部件性能参数变化情况的统计数据[19],选择单部件故障模式:风扇效率变化系数突变-0.04;多部件故障模式:风扇效率变化系数突变-0.04,压气机流量变化系数突变-0.02。仿真时间为10 s,突变故障发生时刻均设定为第2 s,采样步长为20 ms。综合估计精度用全部仿真时间10 s内的均方根误差RMSE来衡量,其表达式为
式中:定义诊断步数为从故障突变到健康参数稳定在故障值±0.01范围内的仿真步数差,用来反映诊断速度;Δh为诊断结果稳定后的平均值,σh为诊断结果的标准差。
在单故障模式下,标准粒子滤波和改进粒子滤波的20次运行平均性能见表2,故障诊断值及其标准差见表3并如图3所示。在多故障模式下,标准粒子滤波和改进粒子滤波的20次运行平均性能见表4,故障诊断值及其标准差见表5并如图4所示。
表2 单故障模式下滤波性能对比
表3 单故障模式下估计值及标准差
图3 单故障模式诊断结果
表4 多故障模式下滤波性能对比
表5 多故障模式下估计值及标准差 %
图4 多故障模式诊断结果
在表3、5中Δh和σh均为20次运行的均值,而σh的均值为8个健康参数的标准差的均值,综合体现诊断结果噪声水平。图3、4中每个健康参数对应的3 个点分别为 Δh-σh,Δh,Δh+σh,图中线段长度反映该健康参数估计值的噪声水平。
由于改进粒子滤波采用更能反映突变情况的伪协方差代替协方差阵,抽样得到的粒子能够更快地逼近真实值,从表2、4中可见,突变故障诊断步数明显减少,2种故障模式下诊断步数均减少约27%;另外由于采用自适应的似然分布函数,提高了先验分布和似然分布的重叠区域,增加了有效粒子数,从表3、5中可见,诊断结果噪声水平显著下降,在2种故障模式下,估计结果的平均标准差均降低约39%。从图3、4中可见,改进粒子滤波器的8个维度的健康参数诊断结果噪声水平均明显下降,估计值更接近真实值。由于改进粒子滤波在故障诊断速度和诊断结果噪声水平这2个方面的改善,综合估计误差明显减小,从表2、4中可见,2种故障模式下RMSE均减小约38%。
标准粒子滤波算法在发动机故障诊断领域有着广泛应用,但该算法存在着突变故障诊断步数多以及诊断结果噪声水平高的问题,本文提出使用伪协方差和自适应似然分布结合的1种改进粒子滤波方法。在发动机故障诊断中的仿真验证表明,相比于标准粒子滤波算法,改进方法能在相同粒子数的前提下,显著减小诊断的均方根误差,并且减少突变故障的诊断步数。综合考虑单故障模式和多故障模式,改进算法能使诊断速度提高约27%,诊断精度提高约38%。优越的性能表明,该改进粒子滤波方法可以进一步推广应用于发动机故障诊断系统中。