高超声速飞行器是飞行速度高于5 Ma的飞行器，这种飞行器飞行时周围的流场通常具有高超声速气体动力学的基本特性，涉及薄激波层、黏性干扰、熵层、高温真实气体效应等现象。因此，高超声速飞行器是复杂的多变量系统，具有强非线性、强耦合、快时变等特性，其黏性效应显著、工作环境条件变化大，容易产生故障，包括执行器故障、传感器故障及结构损伤故障等。因此，高超声速飞行器的故障诊断与调节技术是目前国内外的研究热点和难点^[1]。

故障诊断技术的发展源于建立监控系统的需求，随着设备精密程度以及复杂化程度的提高，对故障诊断的速度和准确性提出了更高的要求。一方面，通过对传统的概率图模型进行改进，能够提高动态故障诊断性能^[2]。另一方面，利用神经网络强大的分类能力对故障模式进行分类与学习，从而实现故障诊断。后者近年来发展迅速，适应性广，适于处理复杂的非线性系统的故障诊断问题。其中，极限学习机(extreme learning machine，ELM)是一种新型算法，其模型是单隐藏层前馈神经网络，与传统的故障诊断方法相比，该方法具有学习速度快及泛化性能好等显著优点^[3]。相比传统的ELM算法，核极限学习机(kernel extreme learning machine, KELM)采用核映射代替随机映射，其实质是通过核函数隐式地将输入空间中低维的线性不可分样本映射到高维甚至无限维的特征空间，使得原空间的非线性可分问题转化为特征空间中的线性可分问题，是机器学习领域中一类非常重要的方法^[4]。

更进一步，将基于知识的故障诊断方法与以遗传算法为代表的优化算法结合，能够有效提高诊断效能^[5]。因此，本文采用KELM结合生物地理学优化(biogeography-based optimization，BBO)算法的方式，开展飞行器故障诊断方法研究。

1 飞行器建模及姿态控制 1.1 飞行器运动学建模

本文的研究对象是美国航空航天署的Langley研究中心公布的Winged-Cone锥形体飞行器模型，经过风洞实验、流体力学仿真和试飞实验数据分析。该型号飞行器获得了相对完备的气动系数表，已有的研究已经将气动参数拟合为非线性函数。同时，公布的资料中也包含了气动、结构等系统的相关参数和公式，是一套成熟的科研模型^[6-7]。飞行器几何形状及其主要的外形尺寸数据如图 1所示。

参考文[8-10]，考虑飞行器的无动力六自由度运动，受重力、空气动力和推力器反作用力作用，非线性微分方程可总结如下：

$ {\dot x = v{\kern 1pt} {\kern 1pt} {\rm{cos}}\gamma {\kern 1pt} {\kern 1pt} {\rm{cos}}\chi ,} $

(1)

$ {\dot y = v{\kern 1pt} {\kern 1pt} {\rm{cos}}\gamma {\kern 1pt} {\kern 1pt} {\rm{sin}}\chi ,} $

(2)

$ {\dot z = - v{\kern 1pt} {\kern 1pt} {\rm{sin}}\chi ,} $

(3)

$ \begin{array}{*{20}{c}} {\dot v = \frac{1}{m}( - D + Y{\kern 1pt} {\kern 1pt} {\rm{sin}}\beta - m{\kern 1pt} {\kern 1pt} {\rm{sin}}\chi + }\\ {{T_{{\rm{r}}x}}{\kern 1pt} {\kern 1pt} {\rm{cos}}\beta {\kern 1pt} {\kern 1pt} {\rm{cos}}\alpha + {T_{{\rm{r}}y}}{\kern 1pt} {\kern 1pt} {\rm{sin}}\beta + {T_{{\rm{r}}z}}{\rm{cos}}\beta {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha ),} \end{array} $

(4)

$ \begin{array}{*{20}{c}} {\dot \gamma = \frac{1}{{mv}}[L{\kern 1pt} {\kern 1pt} {\rm{cos}}\mu - Y{\kern 1pt} {\kern 1pt} {\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{cos}}\beta - m{\rm{cos}}\gamma + }\\ {{T_{{\rm{r}}x}}({\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha + {\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha + {\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\beta {\kern 1pt} {\kern 1pt} {\rm{cos}}\alpha ) - }\\ {{T_{{\rm{r}}z}}({\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{cos}}\alpha - {\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\beta {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha ) - {T_{{\rm{r}}y}}{\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{cos}}\beta ],} \end{array} $

(5)

$ \begin{array}{l} \begin{array}{*{20}{l}} {\dot \chi = \frac{1}{{mv{\kern 1pt} {\kern 1pt} {\rm{cos}}\gamma }}[L{\kern 1pt} {\kern 1pt} {\rm{sin}}\mu + Y{\kern 1pt} {\kern 1pt} {\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{cos}}\beta + }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {T_{{\rm{r}}x}}({\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha - {\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\beta {\kern 1pt} {\kern 1pt} {\rm{cos}}\alpha ) + }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {T_{{\rm{r}}y}}{\kern 1pt} {\kern 1pt} {\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{cos}}\beta - {T_{rz}}({\rm{sin}}\mu {\kern 1pt} {\kern 1pt} {\rm{cos}}\alpha + } \end{array}\\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{sin}}\beta {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha )], \end{array} $

(6)

$ \begin{array}{*{20}{l}} {\dot \alpha = - \frac{1}{{mv{\kern 1pt} {\kern 1pt} {\kern 1pt} {\rm{cos}}\beta }}(L + m{\kern 1pt} {\kern 1pt} {\rm{cos}}\gamma {\kern 1pt} {\kern 1pt} {\rm{cos}}\mu - {T_{{\rm{r}}x}}{\kern 1pt} {\rm{sin}}\alpha + }\\ {{T_{{\rm{r}}z}}{\kern 1pt} {\rm{cos}}\alpha ) + q - {\rm{tan}}\beta (p{\kern 1pt} {\rm{cos}}\alpha + r{\kern 1pt} {\rm{sin}}\alpha ),} \end{array} $

(7)

$ \begin{array}{*{20}{c}} {\dot \beta = p{\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha + r{\kern 1pt} {\rm{cos}}\alpha + \frac{1}{{mv}}(Y{\kern 1pt} {\kern 1pt} {\rm{cos}}\beta + }\\ {m{\kern 1pt} {\kern 1pt} {\rm{cos}}\gamma {\kern 1pt} {\kern 1pt} {\rm{sin}}\mu - {T_{{\rm{r}}x}}{\kern 1pt} {\rm{sin}}\beta {\kern 1pt} {\kern 1pt} {\rm{cos}}\alpha + }\\ {{T_{{\rm{r}}y}}{\kern 1pt} {\rm{cos}}\beta - {T_{{\rm{r}}y}}{\kern 1pt} {\rm{sin}}\beta {\kern 1pt} {\kern 1pt} {\rm{sin}}\alpha ),} \end{array} $

(8)

$ \begin{array}{*{20}{l}} {\dot \mu = p\frac{{{\rm{cos}}\alpha }}{{{\rm{cos}}\beta }} - r\frac{{{\rm{sin}}\alpha }}{{{\rm{cos}}\beta }} + \frac{1}{{mv}}[L({\rm{sin}}\gamma {\rm{sin}}\mu + }\\ {{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\rm{tan}}\beta ) + Y{\kern 1pt} {\kern 1pt} {\rm{sin}}\gamma {\kern 1pt} {\kern 1pt} {\rm{cos}}\mu + m{\kern 1pt} {\kern 1pt} {\kern 1pt} {\rm{cos}}\gamma {\kern 1pt} {\kern 1pt} {\rm{cos}}\mu {\kern 1pt} {\kern 1pt} {\rm{tan}}\beta ],} \end{array} $

(9)

$ {\dot p = \frac{{{I_y} - {I_z}}}{{{I_x}}}rq + \frac{1}{{{I_x}}}({l_{\rm{A}}} + {l_{T{\rm{r}}}}),} $

(10)

$ {\dot r = \frac{{{I_z} - {I_x}}}{{{I_y}}}pq + \frac{1}{{{I_y}}}({m_{\rm{A}}} + {m_{T{\rm{r}}}}),} $

(11)

$ {\dot q = \frac{{{I_x} - {I_y}}}{{{I_z}}}rp + \frac{1}{{{I_z}}}({n_{\rm{A}}} + {n_{T{\rm{r}}}}).} $

(12)

上述各式确定了飞行器在平面大地条件下的质心位置(x、y、z)，单位m；速度v，单位m/s；航迹方位角χ，单位rad；航迹倾角γ，单位rad；攻角α，单位rad；侧滑角β，单位rad；速度滚转角μ，单位rad；三轴角速度(p、q、r)，单位rad/s；一共12个状态的六自由度非线性方程组，同时也包含了气动舵面所提供的气动力(阻力D、升力L、侧向力Y), 单位N；气动力矩(l_A、m_A、n_A)，单位N·m；RCS系统所提供的反作用力之间的关系(T_rx、T_ry、T_rz)，单位N；反作用力矩(l_Tr、m_Tr、n_Tr)，单位N·m；三轴转动惯量(I_x、I_y、I_z)，单位kg·m²；飞行器的质量m，单位kg。

1.2 RCS姿态控制系统 1.2.1 反作用推力器安装布局

高超声速飞行器的飞行过程可以分为起飞段、上升段、入轨段、在轨飞行段、离轨段和无动力返回段等阶段。其中，无动力返回段又可以分为初期再入段、末端区域能量管理段和进场着陆段。从控制角度看，再入段初期是研究的重点，其初始高度为20~70 km，速度为10 Ma~25 Ma。再入段初期大气密度低，导致高超声速飞行器的动压很低，气动舵面的操纵效率不足以满足控制需要，因此通常启动反作用控制系统(reaction control system, RCS)系统协助气动舵共同完成姿态控制任务。

RCS系统采用八喷管推力器模型设定如下：整体呈8个方向均匀分布，编号1、2、3、4号喷管位置分别对应S、W、N、E 4个方向，编号5、6、7、8号喷管位置分别对应WS、WN、EN、ES 4个方向，5、7提供+X方向旋转力矩，6、8提供-X方向旋转力矩，1提供-Y方向旋转力矩，3提供+Y方向旋转力矩，2提供-Z方向旋转力矩，4提供+Z方向旋转力矩。推力器截面到重心长度为10 m，推力设定恒定为10⁷ N。1、2、3、4号喷管提供力矩大小为M₁=-M_y=Fl，M₂=-M_z=Fl，M₃=+M_y=Fl，M₄=+M_z=Fl, 正负号只表示方向。5、6、7、8号喷管安装方式相同，能提供力矩大小为M_x=Fd，M_y=Flsinα，M_z=Flcosα。反作用推力器安装方式如图 2和3所示。

根据以上数据计算得到各喷管工作时对飞行器产生的控制力矩，如表 1所示。

从表 1可知，1号和3号喷管只产生偏航控制力矩，2号和4号喷管只产生俯仰控制力矩，4、6、7、8号喷管可同时产生3个通道的控制力矩。

1.2.2 反作用推力器力矩分配

反作用推力器力矩分配规则采用分轴选择逻辑，分轴选择是一种最直观的喷管选择逻辑，基本思想是：将8个RCS喷管分别分配到俯仰、偏航和滚转3个轴上，特定的喷管只在特定的轴上使用，3轴之间不存在喷管的耦合。这种分配的结果是：每个轴只有一种喷管组合，因此只能产生一种固定大小的力矩。基于上述思想，采用分轴喷管选择方案如下：1号和3号喷管控制偏航轴，2号和4号喷管控制俯仰轴，5、6、7、8号喷管控制滚转轴。具体喷管组合方式及各轴的RCS力矩如表 2所示。

1.2.3 RCS控制下的姿态响应

根据已有的高超声速飞行器六自由度RCS姿态控制系统进行仿真，仿真条件为：H=33.5 km，v= 15 Ma，α=0°，β=0°，μ=0°，p=q=r=0。再入段飞行器需要大迎角以降低飞行速度，故取姿态指令为：α=45°，β=0°，μ=0°。通过仿真得到三轴指令力矩以及根据上述分轴选择规则分配到各喷管的指令力矩。然后，将各喷管指令力矩经过脉宽调制(pulse-width modulation, PWM)后输入到各喷管。力矩分配后的反推力作用器合成力矩接入控制回路，仿真后所得飞行器三轴姿态角如图 4所示，图 5为未分配力矩时的姿态响应。

比较图 4和5可知，使用RCS系统进行姿态控制时，得到的攻角响应基本实现输入指令要求，但是超调量大、响应速度较慢，且达到稳态时干扰较大。总体来说，该RCS姿态控制系统能够通过各推力器的开关，实现飞行器俯仰通道上的姿态控制。

2 RCS推力器故障诊断 2.1 RCS推力器故障建模

由图 4可知，使用该力矩分配方法后飞行器姿态在1.5 s时已经基本达到稳态，为建立飞行器的RCS推力器故障模型，在1.5 s时向其1号喷管注入多种故障：延迟增加、完全失效、推力降低、开关失灵等。加入各故障模型后得到的仿真输出结果如下。

1) 正常工作F′₁=F₁。

正常工作下攻角响应如图 6所示，飞行器攻角在1.5 s之前达到稳态。

2) 延迟增加Δt=1 ms→10 ms (t≥1.5 s)。

延迟增加故障下攻角响应如图 7所示，在1.5 s之后，推力器的延迟增加故障使得飞行器攻角的稳态响应产生了微小的相移。

3) 完全失效F′₁=0 (t≥1.5 s)。

完全失效故障下攻角响应如下。由图 8可知，在1.5 s之后，推力器的完全失效故障使得飞行器攻角状态量失稳，不断增大。

4) 推力降低F′₁=0.2×F₁(t≥1.5 s)。

推力降低故障下攻角响应如图 9所示，在1.5 s之后，推力器的推力降低故障使得飞行器攻角的稳态振荡的周期加长、振幅变大。

5) 开关失灵F′₁=10⁷ N (t≥1.5 s)。

开关失灵故障攻角响应如图 10所示，在1.5 s之后，推力器的开关失灵故障使得攻角状态量失稳，不断减小。

综上所述，除推力器延迟增加故障以及推力降低故障，在其他2种故障的影响下飞行器攻角变化明显。飞行器攻角观测值的变化随着故障类型的不同而不同，故可由攻角数据入手进行基于核极限学习机的RCS推力器故障诊断。

2.2 故障信号特征提取 2.2.1 小波包分解

小波变换具有时频特征的分析功能，是一种线性变换，用频率不同的振荡函数作为窗口函数对信号进行扫频和平移。物理意义上，小波变换系数的模反映的是信号随时间和频率的变化情况。以小波包分解算法为基础的特征提取方法，可以对一定宽度时间窗内的传感器故障数据在频域实现分解。近空间高速飞行器攻角测量信号的特征提取技术，即将传感器输出的原始信号降维，便于进行后续的故障诊断工作。一个时序信号的低频部分是描述信号的主体成分，重构后近似于原始信号，相当于原始信号的消噪；而其高频部分体现出了信号的细节部分，包含了噪声与故障引起的突变成分。因此各部分的能量分布可以有效反映出信号中的突变故障信息^[11-12]。在攻角测量信号中，除了有利用价值的故障信息外，还包含有噪声等与故障无关的信息。故障可能无法通过时域波形判断，利用小波包分解算法处理传感器数据，经过加工提取不同频率段的能量特征信号。

2.2.2 特征提取过程

为了验证基于小波包分解结合KELM故障诊断算法的有效性，在推力器每种故障模式下以不同的飞行状态获得50组数据，包括正常状态共得到250组数据。取RCS模型的0.650~1.650 s时间段，对各组数据进行小波包分解，小波包分解层数设为3层，数据划分成8个频段。对每组数据进行小波包变换，并计算能量特征值。表 3展示了部分训练数据的特征向量，其中第一行数据为故障分类标签(编号由0至4分别对应故障类别为：无故障、延迟增加、完全失效、推力降低、开关失灵)。

从表 3可以看出，每种故障的特征向量有一定的相似性，但通过直观分析数据，难以有效识别故障。任意选取每种故障下的2组数据，可以对比不同故障下的能量特征如图 11和12所示(节点3—8的高频能量值)。

从图 11中可以看出，不同故障的能量特征值在各频段都不尽相同。其中，无故障、延迟增加和推力降低故障的能量特征在高频节点较为接近，但在低频节点相差较大。同时，对比图 11和12可以看出，同一种故障的特征向量在不同组飞行数据中有一定的相似性。

2.3 基于核极限学习机的故障诊断 2.3.1 核极限学习机原理

在隐藏层的特征映射h(x)具体形式未知的情况下，需要引入核函数来度量样本之间的相似度，可以根据Mercer条件定义KELM的核矩阵，KELM的模型的输出表示如下：

$ \begin{array}{l} f(x) = h(x)\beta = h(x){\mathit{\boldsymbol{H}}^ + }{(I/C + \mathit{\boldsymbol{H}}{\mathit{\boldsymbol{H}}^ + })^{ - 1}}\mathit{\boldsymbol{T}} = \\ {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\left[ {\begin{array}{*{20}{c}} {(K(x,{x_1})}\\ \cdots \\ {K(x,{x_N})} \end{array}} \right]^{\rm{T}}}{(I/C + \mathit{\boldsymbol{H}}{\mathit{\boldsymbol{H}}^ + })^{ - 1}}\mathit{\boldsymbol{T}}. \end{array} $

(13)

其中T表示训练样本期望输出矩阵。

从式(13)可以看出，核矩阵Ω_ELM=HH^T仅与输入数据x_i以及训练样本的个数有关。在KELM中，通过核函数K(x_i, x_j)将低维输入空间中的数据(x_i, x_j)转化为高维特征空间中的内积h(x_i)·h(x_j)，与特征空间的维数无关，可以有效避免维数灾难问题。

2.3.2 故障诊断实现过程

1)仿真获取近空间高速飞行器攻角数据。

(1) 建立Winged-Cone模型，加入控制使系统稳定，在攻角指令下仿真飞行，以加入方差为0.000 1°噪声的状态变量作为近空间飞行器攻角传感器的数据。

(2) 分别在t=1.5 s，向一号喷管仿真输出信号注入延迟增加故障、完全失效故障、推力降低故障和开关失灵故障，得到故障发生后的飞行器攻角传感器数据。

2) 极限学习机的训练。

(1) 将上一步中得到的传感器数据经小波包分解后，进一步转化成特征向量，然后平均分成2部分，以其中一部分作为训练样本。

(2) 利用特征向量与核函数，直接获取输出层权重v。矩阵v即是经过训练后得到的映射矩阵，反映了隐藏层到故障的对应关系。

(3) 利用训练样本实现极限学习机的训练。

3) 测试过程。

(1) 以另一半特征向量作为测试样本，输入到训练好的极限学习机中。

(2) 通过极限学习机的输出值，获取故障诊断结果。

2.3.3 仿真结果与分析

选用径向基函数作为核函数，初步设置正则化参数C=1，核参数s=0.1。诊断结果如表 4所示。(故障类型：无故障、延迟增加、完全失效、推力降低和开关失灵)。

表 4结果可以总结为，经过训练和测试过程，训练精度达到93.6%，部分故障类型测试精度达到100%，测试时间均低于0.2 ms。无故障和延迟增加故障的诊断精度较低，其原因可能是延迟增加故障仅使得攻角信号产生了微小的相移，而未对信号的频率特性产生显著影响，故经过小波包分解后的特征数据相对于无故障情况变化不大，相较而言难以诊断。

3 故障诊断方法优化 3.1 核极限学习机参数的优化 3.1.1 生物地理学优化算法

BBO算法建立在生物地理学物种迁移理论上物种迁移和突变而寻求最优解。差分进化(differential evolution, DE)算法也是一种启发式搜索算法，与其他进化算法相同，都包含变异、交叉和选择3个过程。在变异过程中，利用2个个体的差分的加权，与其他个体生成新的个体，然后再进行交叉和选择过程。DE算法是近年来的研究热点，在众多启发式算法中表现出了优异的性能。相比应用广泛的遗传算法和粒子群算法，在计算速度、局部搜索能力等重要性能上，都表现出了更好的效果^[13]。

3.1.2 基于DEBBO的核极限学习机参数优化

径向基函数作为核函数的极限学习机包含2个待优化参数，正则化参数C和核参数s。其中正则化参数C影响训练精度和稳定性，可以认为是平衡训练精度与测试精度的参数。C越大，系统训练精度越高，然而稳定性较差，其取值范围通常选定在[10^-4, 10⁴]。另一个核参数s用于控制核函数的偏差和方差，s越大，则偏差越大，方差越小，取值范围通常选定在[0.1, 10]。

极限学习机的训练结果体现在4个方面：训练精度A_train、测试精度A_test、训练时间和测试时间。参数C和s的选取影响其中的A_train和A_test，因此将优化的目标函数设置为

$ {\rm{min}}F = \mathit{\boldsymbol{K}} \cdot {[{A_{{\rm{ train }}}},{A_{{\rm{ test }}}}]^{\rm{T}}}. $

(14)

其中K=[k₁, k₂]是常值系数矩阵。

使用DEBBO算法优化核极限学习机，优化后的正则化参数C和核参数s为[2 129.7, 6]，使用优化得到的参数训练核极限学习机，进行故障诊断，结果如表 5所示。

比较表 4和5可知，使用DEBBO算法优化核极限学习机参数之后，训练精度和测试精度均有提升至100%，优化效果显著，但延迟增加故障的诊断精度为80%，仍相对较低。

3.2 核函数的优化 3.2.1 单一核函数的选择与比较

常用核函数^[14]如下：

1) 线性核函数

$ K(\mathit{\boldsymbol{x}},\mathit{\boldsymbol{y}}) = {\mathit{\boldsymbol{x}}^{\rm{T}}}\mathit{\boldsymbol{y}}. $

(15)

其中：x、y分别表示故障信号的2组特征向量。

2) 多项式核函数

$ K(\mathit{\boldsymbol{x}},\mathit{\boldsymbol{y}}) = {(\alpha {\mathit{\boldsymbol{x}}^{\rm{T}}}\mathit{\boldsymbol{y}} + c)^d}. $

(16)

其中：d为多项式核函数的参数，不同的d代表不同的函数属性，用来设置多项式核函数的最高项次数。

3) 径向基核函数

$ K(\mathit{\boldsymbol{x}},\mathit{\boldsymbol{y}}) = {\rm{exp}}\left( { - \frac{{{{\left\| {{\kern 1pt} \mathit{\boldsymbol{x}} - \mathit{\boldsymbol{y}}{\kern 1pt} } \right\|}^2}}}{{2{s^2}}}} \right). $

(17)

其中：s为径向基核函数的宽度参数，不同的s代表不同的函数属性。

4) 二层神经网络核函数

$ K(\mathit{\boldsymbol{x}},\mathit{\boldsymbol{y}}) = {\rm{tanh}}(\alpha {\mathit{\boldsymbol{x}}^{\rm{T}}}\mathit{\boldsymbol{y}} + c). $

(18)

其中：tanh为双曲正切函数；α、c分别为核函数比例参数和自由参数，用以控制核函数的空间映射能力。

接下来对比几种核函数在该诊断方法中的使用效果，如表 6—8所示。

通过比较可知，选取线性基核函数和多项式核函数时该方法无法分辨无故障和延迟增加这2种情况，选取sigmoid核函数时该诊断方法完全失效，故单一核函数条件下选取径向基函数最合适。

3.2.2 复合核函数

当目标诊断数据存在异构数据集和多数据源的问题时，只选用一个核函数往往达不到预期效果，若选用多核函数，则可以提高故障诊断的准确性和稳定性^[15]。因此，构造一个多核核函数对提高该方法的智能诊断效果具有非常重要的作用。在满足Mercer定理的前提下，复合核函数通常有3种构造原理^[16]：线性组合、乘积组合、线性组合和乘积组合的叠加。

选择多项式核、径向基核作为复合函数的基函数。根据复合核函数的构造原理以及Mercer准则，构造复合核函数，表示如下^[17]:

$ K(\mathit{\boldsymbol{x}},\mathit{\boldsymbol{y}}) = \frac{1}{2}(\rho {K_{\rm{p}}} + (1 - \rho ){K_{\rm{r}}} + {K_{\rm{p}}}{K_{\rm{r}}}). $

(19)

其中：ρ为复合核函数修正系数，K_p为多项式核函数，K_r为径向基核函数。使用上述复合核函数进行故障诊断的结果如表 9所示。

比较表 5和9可知，每种故障类型的训练精度均有提升，无故障类型测试精度提升了8%，延迟增加故障类型测试精度提升了4%，可见该优化方法同时适用于这几种故障类型。同时，由于复合函数的形式相对复杂，训练时间稍有增长，但测试时间却有缩减，总体来说优化效果显著。

4 结论

本文进行了高超声速飞行器RCS系统推力器故障下的故障在线诊断技术研究，采用了核极限学习机。从故障诊断结果可以看出，应用核极限学习机进行故障诊断可以达到极高的故障诊断精度和较短的测试时间，从理论和仿真的角度验证了基于核极限学习机的飞行器故障诊断方法的可行性。