民机极限飞行状态的动态气动力试验与建模

岑飞，李清，刘志涛，蒋永，张磊

1. 清华大学 自动化系，北京 100084

2. 中国空气动力研究与发展中心 低速空气动力研究所，绵阳 621000

飞行安全是大型商用运输机设计研究中备受关注的核心问题。按照国际民航组织(ICAO)和商业航空安全委员会(CAST)分类标准[1]，飞行事故按类型分为：可控飞行撞地(CFIT)、空中失火(F-NI)、燃油泄漏(FUEL)、飞行失控(Loss of Control-in Flight, LOC-I)、空域冲突(MAC)、冲出跑道(RE)、动力系统故障(SCF-PP)等。按事故类型统计商用航空运输中人员伤亡飞行事故显示[2]，飞行失控在2008—2017的10年间共造成14起事故、1 131人 死亡，在所有飞行事故类型中，无论事故数量还是死亡人数均占比最高。对飞行失控事故的飞行记录数据分析表明[3]，因恶劣天气、系统故障或机组人为因素等使飞机进入超出正常飞行包线的极限飞行状态，是造成飞行失控的重要原因之一。近年发生的波音737MAX事故进一步说明，在一些特定的民用飞机设计约束之下，飞机迎角的快速增加，并进而进入失速或过失速状态，是民用飞机在节能、减排等多种约束目标优化的情况下可能会出现的问题，而目前波音737MAX的处理方法存在很大的风险，因此需要从飞机本体特性和新飞行控制律的角度加以研究。事实上，按线性、定常条件设计的飞机遭遇非线性、非定常的动态气动力环境, 一直存在潜在的风险，如果对其动态气动特性缺乏足够的认识，没有建立精确的动态气动力模型以支撑驾驶员进行应对此种意外情况的应急培训，其结果往往是灾难性的。通过开展民机极限飞行状态的动态气动力特性研究, 来改善飞机飞行失控预防、极限状态改出、飞行模拟训练和飞行事故分析等,以减少或杜绝类似事件再次发生,无疑是今后减少和避免国内外诸多空难事故需要探索的一个重要技术途径。

飞机不慎进入超出正常飞机包线的极限飞行状态时，有两个显著特征[4]：一是迎角、侧滑角变化范围广；二是飞机处于运动快速变化状态。这种复杂的非常规运动环境与人们熟知的定常、直线飞行有着本质区别，这不仅体现在飞机的气动力特性本身——气动力具有完全的动态特征，即非线性、非定常性质；而且体现在气动力风洞试验设计与动态气动特性研究的关系方面，即这种情况下的动态气动特性研究究竟需要怎样的风洞试验数据，以及风洞试验数据如何用于构建气动力模型的问题[5]。近30年来，关于飞机的大迎角与动态气动特性研究，主要围绕战斗机大迎角过失速机动的气动力试验与建模开展，发展起来的动导数、大振幅振荡、旋转天平等动态气动力风洞试验方法[6-7]，已经在战斗机过失速机动与尾旋研究中得到广泛应用。相应的，发展了多种建模方法，包括以气动导数模型[8]、非线性阶跃响应模型[9]、状态空间模型[10]、微分方程模型[11]等为代表的数学建模方法以及以神经网络模型[12]、模糊逻辑模型[13]和支持向量机模型[14]等为代表的人工智能建模方法，这些方法应用在战斗机布局飞机气动特性研究与建模中，使得战斗机大迎角机动过程中的动态气动力问题得到相当充分的认识和不同程度的解决[15-17]。但是，动态气动力特性与飞机布局形式密切关联[18]，而目前对于运输类布局飞机极限飞行状态气动力试验与建模研究仍然有限，对于民机失速，尤其是过失速动态气动特性尚缺乏深入的风洞试验研究；目前飞行模拟器中对于超出正常迎角/侧滑角包线范围的数据，一般是在正常包线数据基础上外推[19]，存在定性的误差，也难以支撑应对极限飞行状态的飞行模拟训练。

本文针对典型民机布局飞机模型，主要关注极限飞行状态中动态气动力的非线性与非定常特征，从风洞动态试验方法中蕴含的假设和应用条件出发，结合极限飞行状态特点和动导数、大振幅试验结果，分析大迎角区域动态气动力参数影响规律和非线性模型结构，在此基础上探讨极限飞行状态动态气动力的风洞试验设计和非定常气动力建模问题。

1 飞机模型

考虑到布局典型性，选择NASA通用运输机标模(Common Research Model, CRM)作为研究模型。CRM是NASA发布的代表典型双发、远程、双通道宽体商用运输机布局的标模，飞机三维数模及数据面向国际合作公开发布[20]。该标模对机翼进行全新设计，采用现代先进的超临界翼型，飞机机翼展弦比为9.0，根稍比为0.275，机翼1/4弦线后掠角为35°；而飞机的机身、平尾和垂尾等部位的关键尺寸、布局参数与波音777-200保持一致[21]，设计巡航马赫数Ma=0.85(设计点升力系数CL=0.5)，飞机气动布局如图1所示。值得一提的是，在NASA发布的CRM原始标模中，飞机不带操纵面，本研究中为了后续开展极限飞行状态操控特性研究需要，参考波音777-200操纵面设计了升降舵、副翼和方向舵。

图1 CRM布局

采用2.45%缩比模型开展动态气动力风洞试验，模型主要参数如表1所示[22]。

表1 CRM动态试验模型参数[22]

2 试验方法

2.1 动态气动力数学模型及其简化形式

飞机受到的气动力可以严格地表示为飞行状态变量的泛函：

Ci(t)=Ci(M(t),h(t),δ(t),α(τ),β(τ),

p(τ),q(τ),r(τ)) -∞≤τ≤t

(1)

式中：i=L,D,Y,Mx,My,Mz，分别表示升力、阻力、侧向力和三轴气动力矩；Ci(t)表示在给定t时刻的气动力和力矩，其不仅与t时刻的飞行状态参数有关，而且与飞行状态参数变化的历史有关。这个表达式能够反映空气流动产生的气动力的非定常、非线性实质，但显然并没有简单的途径可以得到这种泛函的数值[23]。为此，目前的工程实践中通常引入3个假设对这一泛函做如下准定常化、线性化处理[7,23]：

假设1假设较远历史状态对t时刻的气动力影响可以忽略，以飞行状态参数α(τ)为例，在t时刻附近展开为泰勒级数：

(2)

则可以用t时刻的α及其各阶导数值代替式(1)中的α(τ)，其他飞行状态参数作类似处理:

(3)

引入该假设后，气动力变成关于飞行状态参数的准定常、非线性模型，大部分实际问题中，只保留到飞行姿态角速度，就可以保证足够的精度，即

(4)

不失一般性，式(4)可以分解为两部分，即静态气动力Ci1和动态气动力Ci2，Ci1是关于飞机平动运动参数的非线性函数(与飞机姿态角相关)，Ci2是关于飞机转动参数的非线性函数(与飞机姿态角速度相关)：

Ci(t)=Ci1(M,h,δ,α,β)+

(5)

假设2假设动态气动力中，各运动参数对气动力影响可以相互解耦，即

Ci2(M,h,δ,α,β;p)+…

(6)

假设3假设动态气动力中，气动力随飞机姿态角速度是线性变化的，即

(7)

为无因次角速率。

由此，可得

(8)

综上，准定常条件下气动力线性化模型表达式为

(9)

基于上述假设的数学模型以及发展的与之相应的风洞动态试验方法，一般都严格限制在小到中等迎角时气动力非线性、非定常现象很弱的范围，对于极限飞行状态，上述假设是否仍然适用及其所带来的影响，以及如何进行试验设计和建模需要结合特定飞机布局特点进行分析。

2.2 试验内容与方法

结合2.1节分析，尽管在大迎角条件下，上述准定常化、线性化方法，难以严格描述复杂非线性特征，但是实际上任何的非线性气动力现象都是从线性气动力出现异常发展而来的[23]。基于这个认识，针对极限飞行状态动态气动力试验与建模，本文基本思路是：利用已经发展成熟的动导数试验方法，开展包含大迎角在内的不同频率、不同振幅的动导数试验，分析不同迎角区域参数影响规律和动态气动力的非线性特征，建立非线性动态气动力模型；利用大振幅试验结果，对比线性和非线性气动力模型，对动态气动力的非线性特征和非线性建模结果进行分析和验证，并进一步研究气动力的非定常特征，进行非定常气动力建模。

动导数试验通过强迫振荡方式进行，分为俯仰强迫振荡、偏航强迫振荡、滚转强迫振荡试验。试验在中国空气动力研究与发展中心(CARDC)FL-14风洞开展(图2)。试验时，模型绕相应体轴在不同迎角作不同振幅、频率的正弦振荡，按照飞机的短周期或荷兰滚运动模态特征，或者飞机飞行运动中涉及的飞行状态参数典型范围，确定参数模拟范围。试验风速为30 m/s，试验迎角范围为-10°～60°，振荡频率范围为0.5～1.5 Hz，振幅范围为3°～10°。

图2 俯仰振荡和滚转振荡试验

大振幅振荡试验在同一套试验装置上开展，在不同迎角区域作振幅分别为10°、15°、20°的运动，频率范围在0.25～1.25 Hz之间。

3 动态气动力非线性特性

3.1 非线性特征分析

以通过俯仰振荡获取俯仰动导数为例进行分析，在俯仰振荡中，飞机运动的数学描述为

(10)

(11)

图3给出了-4°≤α≤60°范围内不同振荡频率f、振幅A下俯仰阻尼导数试验结果，试验参数如表2所示。

表2 俯仰动导数试验参数

图3 不同频率、振幅俯仰阻尼导数

在常规飞行条件下，利用缩比模型开展动导数试验时，主要考虑减缩频率相似[23]，以表征全尺寸飞机和模型之间关于刚体运动和绕流流体运动在时域上的同时性，形如：

(12)

图4中分别给出小迎角(α0=4°)、中大迎角(α0=16°)、大迎角(α0=40°)下，以振荡频率f=1 Hz、 振幅A=5° 开展动导数试验时，绘制的俯仰力矩系数Cm关于迎角的迟滞回线，图中不同“截止频率”指的是数据处理时的滤波截止频率。例如，截止频率为1 Hz(与强迫振荡频率相同，即通常说的保留1阶量)，即认为绕流流场变化频率与飞机刚体运动频率一致，迟滞回线为标准椭圆。众所周知，基于小扰动线性化假设的动导数就是由该1阶量计算得到的。可以看出，在常规小迎角范围或者在迎角特别大时，保留到高阶量(如6阶)与只留1阶量相比迟滞回线是比较接近的，即表示此时用动导数就可以比较准确地描述飞机的动态气动力，俯仰力矩随俯仰角速度是接近线性变化的；但是，如图4(b)所示，在该中大迎角范围时，气动力存在明显的高阶量，因此，俯仰力矩随俯仰角速度变化有明显的非线性特征，此时，常规的线化导数不能精确反映动态气动力特性。从对气动力的频谱分析也可以看出，对应的数据含有1～6阶量(7阶以上基本没有)，一方面说明了该区域动态气动力存在显著非线性，另一方面该结果也表明对于气动力迟滞回线的计算，保留到6阶(当强迫振荡频率为1 Hz时，截止频率为6 Hz)是比较合理的。

图4 不同平均迎角下强迫振荡时俯仰力矩随迎角变化的迟滞回线

上述结果说明了在10°<α<35°范围内，式(7)所示的线性化假设不成立，必须对动态气动力的试验方法和数学模型进行改进。为了进一步分析动态气动力与角速率之间的关系，图4(d)给出了振幅和频率均不同，但是在强迫振荡平均迎角处(α0=15°)最大俯仰角速率相同的两条俯仰力矩迟滞回线。可见，尽管强迫振荡振幅和频率不同，但在最大角速率相同的位置，迟滞回线是非常接近的。这也是3.2节中把最大角速率引入到非线性动态气动力建模中的一个依据。

3.2 非线性气动力模型

由于在飞机的振荡运动中，角速率不仅与频率相关，而且振幅相关，式(11)定义的最大无因次角速率作为反映飞机动态运动的特征参数，可以同时体现频率和振幅对飞机动态气动特性的影响。为此，在本研究中引入最大无因次角速率作为阻尼导数的一个影响因子，以俯仰阻尼导数为例，即

(13)

用α0=15°、A=10°的俯仰大振幅振荡为例对上述气动力模型进行分析和验证。俯仰大振幅运动如式(10)所示，根据表2的动导数试验结果进行关于角速率的线性或非线性建模，线性数学模型表达式为

(14)

非线性数学模型表达式为

(15)

图5、图6给出了不同大振幅振荡频率(对应于不同的角速率区间)下，俯仰大振幅振荡试验结果和线性、非线性建模计算结果的对比。为了清楚识别图中气动力曲线特征，先对图中曲线形态、标注及其物理意义进行说明：以图5为例，图中静态试验对应飞机俯仰角速度q=0(°)/s；在大振幅振荡运动中，由于动态气动力贡献，出现典型的“迟滞回线”，箭头表示对应曲线中飞机角运动的方向(箭头指向α增大的方向则q为正，反之为负)；对于俯仰力矩回线，气动力曲线随飞机运动沿着逆时针方向变化表明俯仰阻尼导数为负，飞机是动稳定的，反之则动不稳定，显然若回线出现交叉表明飞机动稳定性的突变，预示飞机可能会出现俯仰极限环振荡、机翼摇滚等非线性飞行力学现象。

图5 f=0.25，0.5 Hz时气动力建模结果和大振幅振荡试验结果对比(α0=15°,A=10°)

图6 f=1，1.25 Hz时气动力建模结果和大振幅振荡试验结果对比(α0=15°,A=10°)

4 动态气动力非定常特性

从第3节大振幅试验结果可以看出，加入非线性影响因素后，可以捕捉到飞机主要的动稳定性特征，但仅考虑非线性特性仍不足以精确定量描述运输机在进入极限飞行状态后可能面临的大迎角、快速运动时的气动特性，特别是在飞行失控运动中角速率较大时，气动力将呈现出较强的时间相关性。因此，在具体分析特定快速运动过程中的动态气动力时，还需考虑非定常特性影响。为了描述气动力的非定常特性，这里采用Goman-Khrabrov状态空间建模方法(G-K模型)[10]，通过引入描述流场状态的内部变量即气流分离点位置，来描述分离流流场的动态发展过程，从而建立气动力响应的非定常模型。该模型已经成功应用在战斗机布局飞机非定常气动力建模中[15,24]。

本文以G-K模型以及Fan和Lutze等[25]对该模型的发展为基础，以升力系数为例展开分析，模型表达式为

(16)

式中：

CLα(x)=a1+b1x+c1x2

对于该非定常气动力模型，表征非定常特性的时间常数τ1、τ2的确定是关键。针对该问题，文献[10]推导了通过动导数来辨识时间常数的关系式，从而可以利用动导数试验数据确定时间常数。本文给出了另一种方法，即利用静态测力和大振幅振荡试验结果确定模型时间常数的方法，如下所述。

1)δ、αm、CL0、a1、b1、c1只与静态气动力有关，因此，利用静态测力试验结果，用最小二乘拟合进行上述参数辨识。

2) 利用大振幅振荡试验结果，对τ1、τ2、a2、b2、c2同样采用最小二乘拟合进行参数辨识。

按照这组参数建立非定常气动力数学模型，与静态测力和大振幅试验结果对比如图7、图8所示。值得注意的是，上述参数只用了f=1 Hz时的大振幅试验数据进行模型中与非定常特性有关的参数辨识，图7给出了辨识的非定常气动力模型计算结果与试验结果对比，图8中给出了使用该同一组参数(即τ1、τ2等参数保持不变)其他大振幅运动频率时的非定常气动力模型预测结果，与试验结果均吻合良好。可见，对该民机布局飞机，使用上述非定常气动力模型可以准确描述气动力的非定常特征。

图7 非定常气动力建模与大振幅试验结果对比

图8 非定常气动力模型预测与大振幅试验结果对比

综上，本研究表明，一方面，采用Goman-Khrabrov状态空间建模方法，可以准确描述典型民机布局飞机极限飞行状态下动态气动力的非定常特性，表明该模型适用于运输机布局飞机的非定常气动力建模，同时该模型基于流动机理进行建模，物理意义清晰，从而揭示了民机极限飞行状态下非定常气动力背后的流动结构演化和流动机理；另一方面，该研究为非定常模型时间常数的辨识提供了一种相对简便的通过大振幅振荡试验数据进行确定的方法。总体而言，利用第3节所述的非线性试验设计与建模方法可以获得飞机的主要动态气动力特征，预示出现非线性气动力的迎角范围，在此基础上，再针对特定的飞机运动过程，进一步开展大振幅或多自由度动态试验，进行非定常建模，可更准确获得运输机特定极限飞行状态的动态气动特性。该流程与方法为开展极限飞行状态动态气动力试验设计与建模研究提供了可行的技术途径。

5 结 论

1) 在飞机超出常规迎角的极限飞行状态下，特别是在飞机开始失速到初始尾旋阶段(如对CRM布局而言，10°<α<35°区域)，常规动态气动力试验的线性化、定常化假设不成立，需要考虑动态气动力的非线性与非定常特征。

2) 在动导数试验设计中，除了减缩频率，将一个振荡周期中的最大无因次角速率作为影响因素，形成关于迎角和角速率的非线性动态气动力模型，可以更好地捕捉飞机关于动稳定性的关键特征。

3) 采用Goman-Khrabrov状态空间建模方法，可以准确描述典型民机布局飞机极限飞行状态下动态气动力的非定常特性，表明该模型适用于运输机布局飞机的非定常建模。同时也揭示了民机极限飞行状态下非定常气动力背后的流动结构演化和流动机理。利用大振幅试验数据进行模型参数辨识，可以得到一种相对简单的通过动态试验确定Goman-Khrabrov模型中时间常数的方法。

4) 对于极限飞行状态的动态气动力研究，尽管动导数作为线性气动力参数已经不能用来描述飞机的运动形态，但仍可以预示非线性气动力的迎角范围。在此基础上，再针对特定的飞机运动过程，有针对性地进行大振幅试验及非定常气动力建模，可以更准确捕捉特定极限飞行状态的动态气动特性。因此，所建立的试验平台、方法、数据处理与建模等为民机极限飞行状态动态气动力风洞试验设计与建模研究提供了一个可行的途径。