基于递推最小二乘法的变体飞行器模型参数在线辨识

刘昊东，张庆振，郭云鹤，茅佳雯

(1. 北京航空航天大学自动化科学与电气工程学院，北京 100083； 2. 上海机电工程研究所，上海 201109)

0 引 言

高马赫数飞行器是指飞行速度超过5Ma的有翼或无翼飞行器，如各类导弹、飞机、炮弹等[1]，具有巨大的军事价值和经济价值，是21世纪航空航天领域的一个重要研究方向[2-3]。相比于传统的飞行器，高马赫数飞行器具有以下4个明显的优势：飞行速度快；突防成功率高，拦截困难；射程远，杀伤力大；稍纵即逝，探测难度大[4-5]。基于以上优势，各国都相继对高马赫数飞行器展开了深入的研究，并且取得了一定的研究成果[6]。变体飞行器是指一类能够依据具体的飞行环境和任务要求实时调整外形结构的飞行器[7]。变体飞行器将外形参数作为可控变量，利用外形参数对气动特性的影响来改变飞行器的性能，使其能够适应更宽的飞行空域和速域，从而能够适应更复杂的飞行任务和环境；同时，它在飞行过程中能够根据飞行环境及任务实时调整外形，以获得最优的气动和操纵性能，从而达到减小能耗的目的。高马赫数变体飞行器则是指具有变形能力的高马赫数飞行器，凭借高马赫数以及变形带来的双重优点，大大提升了其适用范围，在未来军事或商业航天中将具有巨大的应用潜力[8-9]。

目前高马赫数变体飞行器大多采用机体/发动机一体化设计技术，这种设计会导致高马赫数变体飞行器的各个状态量之间存在非线性、强耦合等特点。并且，相比于传统飞行器，高马赫数变体飞行器自身的气动特性更为复杂，这些因素就导致高马赫数变体飞行器在制导、控制等方面存在着巨大的困难。同时，高马赫数变体飞行器相比于传统飞行器在飞行过程中马赫数和飞行高度的变化范围更大，环境因素变化情况更为复杂，这就会导致高马赫数变体飞行器中某些模型参数发生较大变化[10-11]。因此，有必要采用在线辨识的方法，实时地获取高马赫数变体飞行器模型的各个参数，为飞行控制系统的在线更新、飞行能力的在线评估以及飞行器故障检测等提供更加准确的模型[12]。

系统辨识是利用系统的输入输出数据获取系统模型的技术，主要包括试验设计、模型辨识、参数估计和系统验证4部分。目前系统辨识中常用的辨识算法有递推最小二乘法(recursive least squares，RLS)[13]、递推极大似然法[14]和卡尔曼滤波方法[15]等。

高马赫数变体飞行器模型的非线性较强、耦合性较高且导致模型阶次较高且未知参数较多，因此高马赫数变体飞行器模型参数辨识比传统飞行器模型参数辨识更加复杂。对于非线性模型的辨识问题，可以通过辨识非线性模型里的线性子系统来得到整个非线性模型的辨识结果。本文选择对高马赫数变体飞行器的动力学模型进行简化，忽略耦合因素，选取纵向小扰动线性化模型，分析设计辨识激励信号，采用经典的递推最小二乘法对多输入多输出的线性运动模型进行辨识，通过仿真验证了递推最小二乘算法的有效性。

1 高马赫数变体飞行器辨识模型

辨识模型的建立是系统辨识的一个重要环节，本章建立了高马赫数变体飞行器的纵向小扰动线性运动模型，并以此模型为辨识模型，对其进行了分解。

1.1 纵向小扰动线性运动模型

在对高马赫数变体飞行器模型进行参数辨识之前首先需要建立辨识模型。考虑到高马赫数变体飞行器的六自由度运动方程过于复杂，多个参数之间存在耦合，辨识难度较大，本文对高马赫数变体飞行器的运动方程进行简化，将飞行器的运动分解为纵向运动和横侧向运动，忽略耦合项。考虑到飞行器的纵向运动模型中包含升力系数、阻力系数、俯仰力矩系数等表征飞行能力的主要参数，本文采用下反角可以变化的高马赫数变体飞行器纵向运动模型为基础建立辨识模型。

通过对高马赫数变体飞行器纵向受力进行分析，可得到式(1)所示微分方程组。

(1)

式中：m为飞行器质量；g为重力加速度；v为速度；X为阻力；Y为升力；α为攻角；ϑ为俯仰角；θ为航迹角；Jz为绕z轴的转动惯量；ωz为俯仰角速度；Mz为俯仰力矩；x,y为沿x轴和y轴的位移。

基于建立纵向运动方程的假设，忽略二阶及三阶小量，同时采用以下基本假设，将纵向扰动运动方程组线性化：

1) 忽略横侧向参数的基准值及纵向角速度。即假定在未扰动运动中，侧向运动参数β、γ、γc、ψ、ψc、ωx、ωy、z，舵偏角δx、δy，纵向运动参数中的ωz、α都是小量，可以在方程中忽略其乘积以及这些参数与其他微量的乘积。此外还假定侧向气动力及其力矩对纵向参数的偏导数在基准状态下均为零。

2) 干扰只改变纵向运动参数，不改变侧向参数，即侧向运动参数的偏量为零。

利用偏导数的简略表示法，可得航向静不稳定飞行器纵向小扰动线性化模型为

(2)

1.2 辨识模型的分析

式(2)为一个MIMO系统，由4个多输入单输出(multiple input single output，MISO)子系统组成。由于这4个MISO子系统的输出都毫无关联，并且每个子系统的结构相同，通过分析研究其中一个MISO子系统的辨识问题就可以得到整个MIMO系统的辨识结果。一个MISO系统可以分解为两个单输入单输出(single input single output，SISO)系统，并且这两个SISO子系统满足叠加原理，通过分别研究每一个SISO子系统，最终可以得到整个MIMO系统的辨识结果。

2 参数辨识方法

目前，最小二乘法在参数辨识领域中是应用最为广泛的方法之一。最小二乘法的改进算法——递推最小二乘法，相比于传统的最小二乘法，由于状态方程非线性较强，可以得到更加精确的辨识结果，有效地解决线性系统的模型参数辨识问题。对于非线性系统，可以通过辨识非线性系统中的线性子系统来获得整个非线性系统的辨识结果。本文采用时域递推最小二乘法对高马赫数变体飞行器的纵向小扰动线性模型进行参数辨识。

2.1 最小二乘法

(3)

基于残差建立最小二乘指标函数

(4)

(5)

取残差展开式的一阶项得

(6)

联立式(4)和式(6)得

(7)

(8)

联立式(7)和式(8)得

(9)

将式(9)简记为

A(χ-χ0)+B=0

(10)

式中，

(11)

(12)

2.2 递推最小二乘法

由于系统的状态方程和观测方程之间的非线性较强，对式(12)进行一次求解难以得到参数的精确值。因此，利用迭代最小二乘方法进行参数辨识，迭代格式为

(13)

式中，

(14)

3 辨识激励信号设计

在对高马赫数变体飞行器进行系统辨识之前，需要设计适合该系统的辨识激励信号，合适的激励信号可以有效地激发系统的各种运动模态，从而获得足够多的信息，来满足系统辨识的需要。不同的输入激励形式会激发高马赫数变体飞行器不同的运动模态。工程上最常见的激励信号为方波信号，包括3-2-1-1方波、偶极方波、2-1-1方波以及阶跃信号等。在相同的基频条件下，相比于正弦信号，方波激励信号的边缘可以产生更宽的信号频谱，是较为理想的系统辨识激励输入，因此本文采用3-2-1-1方波作为激励信号。

3-2-1-1方波信号的特点是幅值呈正负交替出现，宽度比为3∶2∶1∶1，其激励形式及功率谱分布如图1所示。

图1 3-2-1-1方波输入激励信号及其频谱图Fig.1 3-2-1-1 square wave input excitation signal and its spectrum diagram

由图1可以看出，3-2-1-1方波信号的频带较宽，并且在一定的频带范围内具有较高的能量；3-2-1-1方波信号正负交替的特点使飞行器不会过于偏离原有航迹，并且可以有效地激发飞行器的各种运动模态，是飞行器系统辨识较为理想的输入信号。

4 仿真验证

本章在已建立的高马赫数变体飞行器辨识模型的基础上把状态空间转化成传递函数，将得到的传递函数离散化，最后转化为差分方程，对高马赫数变体飞行器进行参数辨识。

4.1 由状态空间模型构建辨识差分方程

根据式(2)所示的已建立的飞行器纵向小扰动线性运动模型，定义

(15)

(16)

则式(2)可以表示为

(17)

式(17)描述的是多输入多输出系统的模型，本文将其分解为8个单输入单输出的子系统，对8个子系统分别进行最小二乘递推辨识，这种思路便于理解和算法实现。

将式(17)所示的纵向小扰动线性状态模型离散化，可以得到其脉冲传递函数为

(18)

传递函数G(z)可由输出量与输入量之比表示，即

(19)

式中：Y(z)为输出量；U(z)为输入量。

可得系统的差分方程为

b13δe(k-3)+b14δe(k-4)+c11δr(k-1)+c12δr(k-2)+c13δr(k-3)+c14δr(k-4)

(20)

α(k)=-a1α(k-1)-a2α(k-2)-a3α(k-3)-a4α(k-4)+b21δe(k-1)+b22δe(k-2)+

b23δe(k-3)+b24δe(k-4)+c21δr(k-1)+c22δr(k-2)+c23δr(k-3)+c24δr(k-4)

(21)

ϑ(k)=-a1ϑ(k-1)-a2ϑ(k-2)-a3ϑ(k-3)-a4ϑ(k-4)+b31δe(k-1)+b32δe(k-2)+

b33δe(k-3)+b34δe(k-4)+c31δr(k-1)+c32δr(k-2)+c33δr(k-3)+c34δr(k-4)

(22)

ωz(k)=-a1ωz(k-1)-a2ωz(k-2)-a3ωz(k-3)-a4ωz(k-4)+b41δe(k-1)+b42δe(k-2)+

b43δe(k-3)+b44δe(k-4)+c41δr(k-1)+c42δr(k-2)+c43δr(k-3)+c44δr(k-4)

(23)

选取式(20)～式(23)为参数辨识模型，采用递推最小二乘法对其差分方程中的未知参数进行辨识。

4.2 仿真实例

本文以采样频率为100 Hz、飞行高度为25 000 m、飞行速度为2 200 m/s的某变下反角飞行器(下反角可以变化的飞行器，下反角为翼梢小翼与机翼基准面的夹角)为例对其纵向小扰动线性运动模型进行了递推最小二乘参数辨识。本文仿真实验中数据来源为文献[9]。

各参数取值为A11=-0.064 1，A12=0.004 3，A13=-0.004 4，A21=-0.040 4，A22=-0.015 9，A41=1.039 2×10-14，A42=8.806 6×10-16，A44=8.103 1×10-21，B21=0.000 437，B22=0.000 063，B41=0.000 809，B42=0.000 422。

升降舵和下反角输入激励信号均选取3-2-1-1方波信号，如图2所示。

图2 输入激励信号Fig.2 Input excitation signal

高马赫数变体飞行器飞行过程中的噪声干扰大部分来自观测噪声，过程噪声的影响相对较小，因此本文在仿真过程中只考虑输出的观测噪声的影响。基于递推最小二乘辨识算法，对式(20)所示的纵向模型子系统进行参数辨识。

在辨识过程中在输出观测量上加入期望为0、方差为10×10-7数量级的白噪声。

图3 攻角系统辨识过程Fig.3 Attack angle system identification process

在加入噪声后，由图3可以看出，在约4 s内，攻角系统参数辨识结果收敛，仿真结果表明递推最小二乘算法计算效率较高、收敛速度较快，能够满足在线辨识算法的时间要求。同时，由表1可知，辨识结果与参数真值之间误差小于6%，辨识结果达到10×10-7数量级，辨识精度较高，基本满足工程要求。

5 结束语

本文研究了基于递推最小二乘法的高马赫数变体飞行器动力学模型中一些未知参数的在线识别和计算方法。首先，根据高马赫数变体飞行器的六自由度运动模型，建立了适用于飞行品质评估、控制律设计的高马赫数变体飞行器纵向小扰动线性运动模型。将多输入多输出系统分解为多个单输入单输出子系统，并以此作为飞行器在线辨识的基础模型。之后，对高马赫数变体飞行器参数辨识所需要的输入激励信号进行了分析与设计，选取了可以有效激发飞行器各种运动模态并且能有效避免飞行器过于偏离原有飞行航迹的3-2-1-1方波信号。最后，以某高马赫数变体飞行器为例，采用递推最小二乘算法，对飞行器纵向小扰动线性运动模型里的单输入单输出线性子系统的在线辨识进行了仿真验证。仿真的过程和结果表明，本文基于递推最小二乘法对飞行器模型参数进行在线辨识的方法收敛速度快、计算效率高，辨识的平均误差小于6%，辨识精度较高，整体辨识精度满足工程要求。