CN110610027A - 一种基于短时数据的航空发动机解析余度计算方法 - Google Patents

Abstract

本发明提出一种基于短时数据的航空发动机解析余度计算方法，属于发动机控制领域。该方法首先确定当前时刻发动机传感器测量信号矩阵的结构；确定当前时刻待估计的传感器测量信号，向测量信号矩阵中载入除当前时刻待估计传感器测量信号外的其它所有传感器测量信号数据，并生成对应的当前时刻观测矩阵；通过迭代估计当前时刻子空间；计算得到当前时刻待估计传感器测量信号的解析余度。本发明利用发动机短时数据进行计算，不需要提前收集大量运行数据；仅依赖于对发动机一般性质的分析，对不同型号、不同健康状况的发动机均具有较好的效果，有很高的应用价值。

Description

一种基于短时数据的航空发动机解析余度计算方法

技术领域

本发明涉及一种基于短时数据的航空发动机解析余度计算方法，属于发动机控制领域。

背景技术

解析余度技术通过传感器测量信号之间的联系来对传感器测量信号进行估计，从而为相应传感器提供冗余备份，提高整个系统的安全性。

随着发动机数字电子控制技术的发展，机载电子设备具备更高的存储和计算能力，使得解析余度技术在航空发动机上的应用成为可能。相较于硬件余度技术，解析余度具有空间占用小，易形成异构余度以及修改维护方便等优点，借助解析余度，通过故障诊断软件，可以对传感器信号的有效性进行检测，识别出故障传感器并进行信号重构，保证发动机的运行安全，为飞行器提供可靠的动力。在这一过程中，解析余度技术所提供的传感器信号估计值发挥了至关重要的作用，对故障诊断的效果有显著影响。

传统的航空发动机传感器解析余度计算方法包括基于模型的方法和基于数据的方法。基于模型的方法如部件级模型、线性模型等自上世纪80年代起得到了广泛研究，但在实际应用时存在以下不足：

1)模型参数的确定需要准确的设计细节或大量的物理实验，上述条件在设计初期难以满足；

2)形参与物参之间的关系不够明确，当出现设计改进或性能退化时模型参数调整较为困难。

近年来，以神经网络和支持向量机为代表的基于数据的计算方法也得到了较多研究，这些方法的不足之处体现在：

1)需要利用大量数据进行训练，在产品开发前期数据量难以满足；

2)对训练范围外的情形泛化能力不足，适用情形有限。

发明内容

本发明的目的在于克服已有技术的不足之处，提出一种基于短时数据的航空发动机解析余度计算方法。本发明利用发动机短时数据进行计算，不需要提前收集大量运行数据；仅依赖于对发动机一般性质的分析，对不同型号、不同健康状况的发动机均具有较好的效果，有很高的应用价值。

本发明提出一种基于短时数据的航空发动机解析余度计算方法，其特征在于，该方法包括以下步骤：

1)确定当前时刻发动机传感器测量信号矩阵S的结构；定义矩阵S规模为m×n，m为矩阵S中测量信号的采样时刻总数，n代表当前时刻传感器测量信号的总数；

矩阵S的表达式如式(1)所示，矩阵S中每一列代表同一个传感器在不同时刻的测量信号，每一行代表同一时刻采集到的不同传感器的测量信号，其中矩阵第m行代表当前时刻采集的不同传感器的测量信号；

其中，s_i,j表示第j个传感器在第i个时刻测量信号；

2)确定当前时刻待估计的传感器测量信号，向矩阵S中载入除当前时刻待估计传感器测量信号外的其它所有传感器测量信号数据，并生成对应的当前时刻观测矩阵；具体步骤如下：

2-1)确定当前时刻待估计的p个传感器测量信号，1≤p＜n；

2-2)按照步骤1)确定的矩阵S的结构向发动机传感器测量信号矩阵S中载入前m-1个时刻的所有传感器测量信号，以及当前时刻除p个待估计传感器测量信号外的其他传感器测量信号；其中，当前时刻p个待估计传感器测量信号对应位置的元素为未知；

2-3)构建当前时刻观测矩阵H，所述观测矩阵H与测量信号矩阵S规模相同；当前时刻观测矩阵H表达式如下：

其中，h_i,j为矩阵H中第i行第j列的元素，h_i,j＝0表示式(1)中发动机传感器测量信号矩阵S的元素s_i,j的值未知，h_i,j＝1表示式(1)中发动机传感器测量信号矩阵S的元素s_i,j的值已知；

3)估计当前时刻子空间U；具体步骤如下：

3-1)当前时刻子空间U初始化；

令U∈R^m×d，d＝min{n_u+5,n}，其中n_u为n个传感器测量信号中用于测量发动机输入状态的传感器测量信号的数目；

在当前时刻为初始时刻或当前时刻步骤1)确定的矩阵S结构相比上一时刻出现变化时，随机生成一个正交矩阵作为当前时刻T第一步迭代开始时的U的初始值；否则在当前时刻T第一步迭代开始时的U的初始值为上一时刻迭代结束时的子空间U的估计值其中sk_max表示每个时刻的迭代次数上限；

3-2)令迭代步数K＝1；

3-3)在当前时刻第K步迭代时，在步骤2)载入传感器测量信号的矩阵S中任意选取一列，假设选择矩阵S的第k列，记为σ(k)，即：

通过式(3)，得到第k列对应的坐标向量ω(k)：

其中，表示当前时刻T第K步迭代的子空间U的初始值，K＝1时，根据步骤3-1)确定，K>1时，为当前时刻T第K-1步迭代后得到的子空间U的更新值；P_η(k)()是一个投影函数，η(k)为观测矩阵H中第k列向量，即：

对任意向量ξ∈R^m×1，函数P_η(k)()定义为：

其中，ξ_i表示向量ξ的第i个元素；

3-4)根据式(5)，计算第K步迭代的更新步长θ(k)，其中，ρ(k)表示本次迭代计算的预测向量，δ(k)表示预测向量与已知值的残差，δ(k)_i表示向量δ(k)的第i个元素；

3-5)更新第K次迭代的子空间

3-6)令K＝K+1，然后重新返回步骤3-3)，每次随机抽取矩阵S的任一列进行子空间U的迭代更新，直到到达迭代次数的上限sk_max，迭代结束，得到当前时刻子空间U的估计值然后进入步骤4)；

4)对步骤2)确定的每个待估计传感器测量信号进行解析余度计算，得到当前时刻待估计传感器测量信号的解析余度计算结果；具体步骤如下：

4-1)在步骤2)载入传感器测量信号的矩阵S中任意选取一列待估计传感器测量信号所在的列向量，假设选择矩阵S的第g列，记为σ(g)，按照式(7)得到该列向量对应的坐标向量ω(g)：

4-2)按照式(8)对步骤4-1)选取的列向量中待估计传感器测量信号值进行估计，得到该传感器测量信号在当前时刻的解析余度

其中，表示向量中第m个元素；

当前时刻该传感器测量信号的解析余度计算完成；

4-3)重复步骤4-1)至4-2)，直至步骤2)确定的所有当前时刻待估计测量信号的解析余度计算完成；

5)当下一时刻到来时，将该时刻作为新的当前时刻，然后重新返回步骤1)。

本发明的特点及有益效果在于：

1)本发明基于航空发动机的一般性质，对航空发动机的性能与健康状况不敏感，对航空发动机工作点变化和性能波动具有很好的鲁棒性；

2)本发明算法基于短时数据，不需要大量数据训练，同时算法为在线算法，数据易得，具有良好的应用前景；

3)本发明算法不依赖特定的输入输出关系，可以为多个传感器信号进行估计，适合的使用场景更多。

附图说明

图1为本发明方法的总体流程图。

图2为本发明一个实施例的低压转子转速N₁传感器解析余度计算效果图。

图3为本发明一个实施例的高压转子转速N₂传感器解析余度计算效果图。

具体实施方式

本发明提出一种基于短时数据的航空发动机解析余度计算方法，下面结合附图和具体实施例对本发明进一步详细说明如下。

本发明提出一种基于短时数据的航空发动机解析余度计算方法，整体流程如图1所示，该方法包括以下步骤：

1)确定当前时刻发动机传感器测量信号矩阵S的结构，包括矩阵S的规模及矩阵中各个位置元素的物理含义。

定义矩阵S规模为m×n，m为矩阵S中测量信号的采样时刻总数，取值范围为20～50，对应的采样频率为2～20Hz，n代表当前时刻传感器测量信号的总数，由发动机控制器接收到的传感器测量信号数目确定，所述的测量信号包括发动机的转速、各截面温度、各截面压力以及发动机执行机构状态，不包含开关量信号。不同时刻矩阵S对应的传感器数量、排列和测量对象可以不同。

矩阵S的表达式如式(1)所示，矩阵S中每一列代表同一个传感器在不同时刻的测量信号，不同传感器测量信号之间排列顺序任意；每一行代表同一时刻采集到的不同传感器的测量信号，其中矩阵第m行代表当前时刻采集的不同传感器的测量信号，第m-1行则代表上一时刻采集到的传感器的测量信号，以此类推。

其中，s_i,j表示第j个传感器在第i个时刻测量信号。

2)确定当前时刻待估计的传感器测量信号，向矩阵S中载入除当前时刻待估计传感器测量信号外的其它所有传感器测量信号数据，并生成对应的当前时刻观测矩阵。具体步骤如下：

2-1)确定当前时刻待估计的p个传感器测量信号，待估计传感器测量信号指的是步骤1)中n个传感器测量信号中某一个(或多个)，除不能选择步骤1)中全部n个传感器测量信号作为待估计对象以外，待估计对象可以任意选取，即1≤p＜n。

2-2)按照步骤1)确定的矩阵S的结构在发动机传感器测量信号矩阵S中载入前m-1个时刻的所有传感器测量信号，以及当前时刻除p个待估计传感器测量信号外的其他传感器测量信号；其中，当前时刻p个待估计传感器测量信号对应位置的元素为未知，其他位置元素均为已知。

2-3)构建当前时刻观测矩阵H，所述观测矩阵H与测量信号矩阵S规模相同，H中值为1的元素表示S中对应位置的元素值已知，H中值为0的元素表示S中对应位置的元素值未知。当前时刻观测矩阵H表达式如下：

其中，h_i,j为矩阵H中第i行第j列的元素，h_i,j＝0表示式(1)中发动机传感器测量信号矩阵S的元素s_i,j的值未知，h_i,j＝1表示式(1)中发动机传感器测量信号矩阵S的元素s_i,j的值已知。

3)估计当前时刻子空间U；子空间U以矩阵形式表示，其具体步骤如下：

3-1)当前时刻子空间U初始化。

令U∈R^m×d，d＝min{n_u+5,n}，其中n_u为步骤1)所述n个传感器测量信号中用于测量发动机输入状态的传感器测量信号的数目。

在当前时刻为初始时刻或当前时刻步骤1)确定的矩阵S结构相比上一时刻出现变化时时，随机生成一个正交矩阵作为当前时刻T第一步迭代开始时的U的初始值；否则在当前时刻T第一步迭代开始时的U的初始值取上一时刻(T-1时刻)迭代结束时的子空间U的估计值sk_max表示每个时刻的迭代次数上限。

3-2)令迭代步数K＝1；

步骤3-3)-3-5)为循环迭代计算过程，下面以迭代计算的第K步为例，介绍每次迭代计算的过程。

3-3)在当前时刻第K步迭代时，随机选择步骤2)载入传感器测量信号的矩阵S的任一列，不失一般性，假设所选择的列向量为矩阵S的第k列，记为σ(k)，即：

通过式(3)，得到第k列对应的坐标向量ω(k)：

其中，表示当前时刻T第K步迭代的子空间U的初始值，K＝1时，根据步骤3-1)确定，K>1时，为当前时刻T第K-1步迭代后得到的子空间U的更新值；P_η(k)()是一个投影函数，η(k)为观测矩阵H中第k列向量，即：

对任意向量ξ∈R^m×1，函数P_η(k)()被定义为:

其中，ξ_i表示向量ξ的第i个元素。

3-4)根据式(5)，计算第K步迭代的更新步长θ(k)，其中，ρ(k)表示本次迭代计算的预测向量，δ(k)表示预测向量与已知值的残差，δ(k)_i表示向量δ(k)的第i个元素。

3-5)更新第K次迭代的子空间

3-6)令K＝K+1，然后重新返回步骤3-3)，每次迭代随机抽取矩阵S的任一列进行子空间U的迭代更新，到到达迭代次数的上限；迭代次数上限sk_max的参考值为d²log(m)，迭代结束，得到当前时刻子空间U的估计值然后进入步骤4)。

4)对步骤2)确定的每个待估计传感器测量信号进行解析余度计算，得到当前时刻待估计传感器测量信号的解析余度计算结果。具体步骤如下：

4-1)在步骤2)载入传感器测量信号的矩阵S中任意选取一列待估计传感器测量信号所在的列向量，假设本步骤选择矩阵S的第g列，记为σ(g)，按照式(7)得到该列向量对应的坐标向量ω(g)。

4-2)按照式(8)对步骤4-1)选取的列向量中待估计传感器测量信号值进行估计，从而得到该传感器测量信号在当前时刻的估计值(解析余度)该测量信号的解析余度计算完成。其中，表示向量中第m个元素。至此，当前时刻该测量信号的解析余度计算完成。

4-3)重复步骤4-1)至4-2)，直至步骤2)确定的所有当前时刻待估计测量信号的解析余度计算完成。

5)当下一时刻到来时，将该时刻作为新的当前时刻，然后重新返回步骤1)。

下面结合一个具体实施例对本发明进一步详细说明如下。

本实施例提出一种基于短时数据的航空发动机解析余度计算方法，该方法包括以下步骤：

1)确定当前时刻发动机传感器测量信号矩阵S的结构，包括矩阵S的规模及矩阵中各个位置元素的物理含义。

定义矩阵S规模为m×n，m＝30，n＝8。本实施例中矩阵S的形式如式(9)所示，矩阵S中每列元素代表的传感器依次为：低压转子转速N₁、高压转子转速N₂、燃油流量W_f、排气温度EGT、环境总温T_t0、环境总压P_t0、滑油压力OIP、滑油温度OIT，为便于表示，用σ(1),σ(2),…,σ(8)分别代表以上8个传感器的信号；矩阵S共包含30个时刻的传感器信号，采样频率率为8Hz，当前时刻传感器测量信号矩阵S结构为：

2)确定当前时刻待估计的传感器测量信号，向矩阵S中载入除当前时刻待估计传感器测量信号外的其它所有传感器测量信号数据，并生成对应的当前时刻观测矩阵。具体步骤如下：

2-1)确定当前时刻待估计的p个传感器测量信号；

本实施例选择低压转子转速N₁、高压转子转速N₂为待估计传感器测量信号。

2-2)按照步骤1)确定的矩阵S的结构在测量信号矩阵S中载入前m-1个时刻的所有传感器测量信号，以及当前时刻除p个待估计传感器测量信号外的其他传感器测量信号；其中，当前时刻p个待估计传感器测量信号对应位置的元素为未知，其他位置元素均为已知。

2-3)构建一个观测矩阵H，所述观测矩阵H与测量信号矩阵S规模相同，H中值为1的元素表示S中对应位置的元素值已知，H中值为0的元素表示S中对应位置的元素值未知。本实施例中，当前时刻矩阵H为：

3)估计当前时刻子空间U，具体步骤如下：

3-1)初始化当前时刻子空间U。

在当前时刻为初始时刻或当前时刻步骤1)确定的矩阵S结构相比上一时刻出现变化时时，随机生成一个正交矩阵作为当前时刻T第一步迭代开始时的U的初始值；否则在当前时刻T第一步迭代开始时的U的初始值取上一时刻(T-1时刻)迭代结束时的子空间U的估计值sk_max表示每个时刻的迭代次数上限。

3-2)令迭代步数K＝1；

步骤3-3)-3-5)为循环迭代计算过程，下面以迭代计算的第K步为例，介绍每次迭代计算的过程。

3-3)在当前时刻第K步迭代时，随机选择步骤2)载入传感器测量信号的矩阵S的任一列，不失一般性，假设所选择的列向量为矩阵S的第k列，记为σ(k)。通过式(3)，得到对应的坐标向量ω(k)。

3-4)通过式(5)，计算更新步长θ(k)。

3-5)通过式(6)，更新子空间

3-6)令K＝K+1，然后重新返回步骤3-3)，每次迭代时随机抽取步骤2)载入传感器测量信号的矩阵S的任一列进行子空间U的迭代更新，直到到达迭代次数的上限sk_max＝50，迭代结束，得到当前时刻子空间估计值然后进入步骤4)。

4)对步骤2)确定的每个待估计传感器测量信号进行解析余度计算，得到当前时刻待估计传感器测量信号的解析余度计算结果。具体步骤如下：

4-1)在测量信号矩阵S中任意选取一列待估计传感器测量信号所在的列向量，本实施例中首先选择低压转子转速N₁传感器所在的列，即σ(1),，按照式(7)得到对应的坐标向量ω(1)。

4-2)按照式(8)对步骤4-1)选择的低压转子转速N₁传感器测量信号值进行估计，从而得到该传感器测量信号在当前时刻的估计值(解析余度)低压转子转速N₁传感器解析余度计算完成。

4-3)选择测量信号矩阵S中高压转子转速N₂传感器所在的列，即σ(2)，重复步骤4-1)至4-2)，得到该传感器测量信号在当前时刻的估计值(解析余度)高压转子转速N₂传感器解析余度计算完成。

5)当下一时刻到来时，将该时刻作为新的当前时刻，然后重新返回步骤1)。

从而得到一段时间内各个时刻的低压转子转速N₁传感器的解析余度计算结果，如图2所示，以及高压转子转速N₂传感器的解析余度计算结果，如图3所示。

在图2中，横坐标代表采样时刻，单位为125ms，纵坐标代表低压转子转速N₁的值，单位为％，图中黑色连续实线代表了该传感器的实际测量信号，黑色点划线代表了本实施例中该传感器的解析余度计算结果，图2表明本实施例中的低压转子转速N₁传感器的解析余度计算结果较为接近实际值，且对发动机的工作状态变化具有较好的适应能力。

在图3中，横坐标代表采样时刻，单位为125ms，纵坐标代表高压转子转速N₂的值，单位为％，图中黑色连续实线代表了该传感器的实际测量信号，黑色点划线代表了本实施例中该传感器的解析余度计算结果，图3表明本实施例中的高压转子转速N₂传感器的解析余度计算结果较为接近实际值，且对发动机的工作状态变化具有较好的适应能力。

Claims (1)

1.一种基于短时数据的航空发动机解析余度计算方法，其特征在于，该方法包括以下步骤：

1)确定当前时刻发动机传感器测量信号矩阵S的结构；定义矩阵S规模为m×n，m为矩阵S中测量信号的采样时刻总数，n代表当前时刻传感器测量信号的总数；

矩阵S的表达式如式(1)所示，矩阵S中每一列代表同一个传感器在不同时刻的测量信号，每一行代表同一时刻采集到的不同传感器的测量信号，其中矩阵第m行代表当前时刻采集的不同传感器的测量信号；

其中，s_i,j表示第j个传感器在第i个时刻测量信号；

2)确定当前时刻待估计的传感器测量信号，向矩阵S中载入除当前时刻待估计传感器测量信号外的其它所有传感器测量信号数据，并生成对应的当前时刻观测矩阵；具体步骤如下：

2-1)确定当前时刻待估计的p个传感器测量信号，1≤p＜n；

2-2)按照步骤1)确定的矩阵S的结构向发动机传感器测量信号矩阵S中载入前m-1个时刻的所有传感器测量信号，以及当前时刻除p个待估计传感器测量信号外的其他传感器测量信号；其中，当前时刻p个待估计传感器测量信号对应位置的元素为未知；

2-3)构建当前时刻观测矩阵H，所述观测矩阵H与测量信号矩阵S规模相同；当前时刻观测矩阵H表达式如下：

其中，h_i,j为矩阵H中第i行第j列的元素，h_i,j＝0表示式(1)中发动机传感器测量信号矩阵S的元素s_i,j的值未知，h_i,j＝1表示式(1)中发动机传感器测量信号矩阵S的元素s_i,j的值已知；

3)估计当前时刻子空间U；具体步骤如下：

3-1)当前时刻子空间U初始化；

令U∈R^m×d，d＝min{n_u+5,n}，其中n_u为n个传感器测量信号中用于测量发动机输入状态的传感器测量信号的数目；

在当前时刻为初始时刻或当前时刻步骤1)确定的矩阵S结构相比上一时刻出现变化时，随机生成一个正交矩阵作为当前时刻T第一步迭代开始时的U的初始值；否则在当前时刻T第一步迭代开始时的U的初始值为上一时刻迭代结束时的子空间U的估计值其中sk_max表示每个时刻的迭代次数上限；

3-2)令迭代步数K＝1；

3-3)在当前时刻第K步迭代时，在步骤2)载入传感器测量信号的矩阵S中任意选取一列，假设选择矩阵S的第k列，记为σ(k)，即：

通过式(3)，得到第k列对应的坐标向量ω(k)：

其中，表示当前时刻T第K步迭代的子空间U的初始值，K＝1时，根据步骤3-1)确定，K>1时，为当前时刻T第K-1步迭代后得到的子空间U的更新值；P_η(k)()是一个投影函数，η(k)为观测矩阵H中第k列向量，即：

对任意向量ξ∈R^m×1，函数P_η(k)()定义为：

其中，ξ_i表示向量ξ的第i个元素；

3-4)根据式(5)，计算第K步迭代的更新步长θ(k)，其中，ρ(k)表示本次迭代计算的预测向量，δ(k)表示预测向量与已知值的残差，δ(k)_i表示向量δ(k)的第i个元素；

3-5)更新第K次迭代的子空间

3-6)令K＝K+1，然后重新返回步骤3-3)，每次随机抽取矩阵S的任一列进行子空间U的迭代更新，直到到达迭代次数的上限sk_max，迭代结束，得到当前时刻子空间U的估计值然后进入步骤4)；

4)对步骤2)确定的每个待估计传感器测量信号进行解析余度计算，得到当前时刻待估计传感器测量信号的解析余度计算结果；具体步骤如下：

4-1)在步骤2)载入传感器测量信号的矩阵S中任意选取一列待估计传感器测量信号所在的列向量，假设选择矩阵S的第g列，记为σ(g)，按照式(7)得到该列向量对应的坐标向量ω(g)：

4-2)按照式(8)对步骤4-1)选取的列向量中待估计传感器测量信号值进行估计，得到该传感器测量信号在当前时刻的解析余度

其中，表示向量中第m个元素；

当前时刻该传感器测量信号的解析余度计算完成；

4-3)重复步骤4-1)至4-2)，直至步骤2)确定的所有当前时刻待估计测量信号的解析余度计算完成；

5)当下一时刻到来时，将该时刻作为新的当前时刻，然后重新返回步骤1)。