郑无计， 李颖晖， 屈亮， 徐浩军， 袁国强

空军工程大学 航空航天工程学院， 西安 710038

基于正规形法的结冰飞机着陆阶段非线性稳定域

郑无计， 李颖晖*， 屈亮， 徐浩军， 袁国强

空军工程大学 航空航天工程学院， 西安 710038

结冰导致飞机飞行包线缩小、对飞行安全产生严重威胁，研究结冰后飞机非线性稳定域对边界保护系统的设计和飞行安全的提高极其重要。以某型运输机为研究对象，考虑飞机非线性气动特性建立飞机纵向非线性模型并进行增稳控制补偿设计；然后通过流形和正规形理论刻画结冰飞机纵向非线性稳定边界并得到稳定边界的解析表达式，通过仿真的手段验证了正规形理论确定的稳定边界和解析表达式的有效性和准确性。最后，分析了飞机着陆过程中，结冰因子对结冰飞机稳定域的影响以及结冰飞机发生事故的机理。研究结果表明，轻度结冰使飞机非线性稳定域缩小；重度结冰导致飞机稳定性发生改变；在未察觉飞机结冰的情况下，飞行员的常规操纵会使飞行状态超出结冰飞机的非线性稳定域、导致飞行事故。研究结果可为飞机结冰后的边界保护提供一定的参考。

飞机结冰； 稳定域； 非线性； 正规形理论； 运输机

飞机结冰是飞机在结冰条件下飞行时，大气中的液态水在部件表面结冰并累积成冰的一种物理过程，在飞行实践中广泛存在[1]。飞机结冰破坏了飞机气动特性，使升力降低、阻力增加，且使飞机失速迎角大幅度降低，恶化了飞机的飞行性能，对飞行安全产生严重威胁。据美国联邦航空管理局(FAA)和美国国家航空航天局(NASA)统计，在1976—1994年间，发生的飞行事故中有近16起与飞机结冰有关，并导致139人死亡[2]。美国Safety Advisor[3]对1990—2000年的飞行事故进行了详细的统计，结果表明由气象因素引起的飞行事故中仅结冰就占12%。1994年，美国鹰航公司的一架ATR72-212飞机在印第安纳Roselawn地区由于机翼结冰而坠毁；随后2002年在中国台湾坠毁的一架货机以及2006年在安徽失事的一架中国军用运输机都是由于飞机结冰引起的。

国外对结冰现象的研究开始于20世纪30年代。NASA在2000年关于结冰对现代飞机影响[4]的报告中系统地研究了二维翼型在不同结冰情况的气动特性。Miller和Ribeens[5]通过对飞行数据进行气动辨识的方法得到了结冰飞机的气动特性，并初步研究了结冰飞机的运动特性。Bragg等[6]研究了结冰飞机飞行动力学特性。国内以研究飞机结冰及防除冰机理[7]作为重点内容，近年来，相关单位开始对结冰飞机的动力学特性进行研究[8-9]。

对于结冰飞机已经开展的主要工作有结冰预警、防除冰以及结冰后边界保护等，而对于结冰飞机而言为确保飞行安全,飞机结冰后进行边界保护措施是必要的，而结冰条件下飞机的飞行安全边界会缩小[10]，根据实际结冰情况确定边界范围是进行边界保护的重要前提，该边界范围即结冰飞机的非线性稳定域，它表示飞机动态稳定性能，反映了某一飞行状态下的抗扰动性，因此开展结冰飞机稳定域的研究具有重要意义。目前国内外对于结冰飞机动力学特性的研究主要集中在结冰因子[6,11]对飞机静稳定性的影响，而对于结冰飞机非线性稳定域的研究较少，高浩和周志强[12]，何植岱和郭文[13]等利用分支突变理论研究了平衡点变化对飞机运动的影响，但研究内容未涉及到飞机的稳定域问题。Robert[14]和詹浩[15]等利用优化控制和不变集理论研究得到了一定时间内的飞行最大可控边界集，但该种方法使用的初始状态空间为一三维长柱形，难以对飞机飞行范围进行精确地描述。

为得到结冰对飞行特性的影响以及飞机失速的机理，本文首次将正规形理论应用到飞机稳定域相关问题的研究，并给出高阶多项式形式的非线性映射的求解方法，具体研究内容如下：研究基于某型运输机，由于该型飞机采用放宽静稳定度技术，因此本文在考虑气动特性非线性建立飞机非线性动力学模型的同时采用增稳控制补偿模型。然后，基于正规形[16-17](Normal Form, NF)理论，给出便于计算机实现的高阶非线性映射求解系统稳定边界的方法并提出利用边界参数w(x)w(xsep)判断系统状态是否位于稳定域内的相关方法，利用高阶正规形法(后文简写为高阶NF法)求解了系统的稳定边界，并对方法的有效性和正确性进行了验证。在此基础上分析了结冰因子变化对飞机稳定边界的影响，并根据其对稳定域的影响程度定义了飞机轻度结冰和重度结冰的结冰因子范围。最后利用稳定边界及参数w(x)w(xsep)研究分析了着陆阶段拉平减速过程中飞机轻度结冰和重度结冰情况对飞行安全的威胁，得到了在飞行员未察觉飞机结冰的情况下，对积冰飞机进行常规操纵也会导致飞行事故的原因。

1 结冰飞机非线性模型的建立

为简化分析过程并得到较为准确的结果，本文采用如下形式的三维动力学模型。

1.1 飞机动力学模型

飞机结冰后的运动方程为状态向量x、控制变量δe以及结冰程度参数η的函数，本文采用的动力学方程可表示为

(1)

式中:状态向量x包括飞行速度Vt、迎角α及俯仰角速度q,即

x=[Vtαq]T

(2)

其中：δe为升降舵偏角。

飞机三维(Vt、α和q)非线性常微分刚体运动方程可表示为

(3)

式中：θ为俯仰角；D为飞机阻力；L为飞机升力；M为飞机纵向力矩；m为飞机质量；Jy为飞机纵向转动惯量。飞机阻力、升力及力矩可表示为

(4)

1.2 多项式形式的空气动力学模型

本文主要研究运输机低速着陆阶段，且该阶段速度变化较小，因此为简化分析暂不考虑飞行马赫数对飞机气动特性的影响，但为了更好地体现低速着陆阶段的特点，应考虑增升装置对运输机气动特性的影响，故选择襟翼偏角为30°的情况进行仿真模拟。根据飞行数据进行拟合[5]可得到较为准确的气动数据随迎角α、俯仰角速度q以及升降舵偏角δe变化的多项式形式的解析表达式为

(5)

式中：多项式系数xi(i=1,2,3,4)、zi(i=1,2,3,4,5)和mi(i=1,2,3,4,5)为结冰程度参数η的函数，根据文献[6,11]所述，可通过式(6)体现结冰程度参数η对气动特性的影响。

CA,iced=(1+ηficed)CA

(6)

式中：CA和CA,iced分别为结冰前、后飞机气动导数值；ficed为结冰系数，反映CA对结冰的敏感性，对于给定的飞机为常值。结冰影响模型式(6)是由Bragg等[6]通过对大量实验数据进行分析、总结得到的，该模型具有简单的表达形式，可以根据不同飞机的结构及飞行条件较为准确地反映出结冰对飞机气动特性的影响情况，可用于研究结冰对飞机性能及品质的影响，且具有一定的通用性。

本文使用的具体参数如表1所示。根据式(5)可以得到结冰前该型运输机纵向相关气动导数，在本文飞行条件下该型运输机结冰前纵向相关气动导数如表1所示，不同结冰情况的相关气动导数可根据式(6)进行计算。

表1结冰前后相关气动导数

Table1Aerodynamicderivativesforcleanandicedaircraft

AerodynamicderivativeCleancaseficedCLα7．8446-0．10Cmα-1．8661-0．5Cmq-38．1058-0．1754

1.3 增稳控制补偿模型

由于本文研究的背景飞机放宽了静稳定度，所以需要对飞机本体的操纵稳定品质进行增稳控制补偿设计，拟采用状态变量迎角α、俯仰角速度q反馈实现飞机结冰前后飞机系统的增稳控制补偿。

δe=kαΔα+kqΔq

(7)

式中:kα和kq为状态反馈系数；Δα=α-α0，Δq=q-q0为系统状态误差；α0和q0分别为飞机本体平衡状态的迎角和俯仰角。

2 正规形理论对稳定域的估计

正规形理论以流形理论[18-19]估计稳定域方法为基础，通过多项式形式的非线性映射实现复杂非线性系统向简单线性系统的变换，并结合线性系统流形的空间拓扑结构特点给出原非线性系统稳定边界的一种近似估计，且该近似边界具有多项式形式的解析表达式。

2.1 流形理论估计稳定域方法

流形理论是根据微分系统的空间拓扑结构确定系统稳定平衡点(Stable Equilibrium Point，SEP)吸引区的一种方法，它利用直接积分方法得到边界上不稳定平衡点(Unstable Equilibrium Point，UEP)的稳定流形作为稳定边界的一种精确估计。

文献[18]指出，微分动力系统在满足条件1～条件3时，系统在稳定平衡点的稳定边界由边界上不稳定平衡点的稳定流形构成。

条件1边界上所有的平衡点为双曲形式的不稳定平衡点，双曲平衡点即系统在该点线性化后特征值实部不为零。

条件2边界上不稳定平衡点的稳定流形与不稳定流形保持横截条件。

条件3边界上的所有轨线最终趋近于特定平衡点。

下面给出基于流形理论估计系统稳定边界的一般方法。

考虑式(8)所示的非线性系统

(8)

式中：x=[x1x2…xN]T。

1) 通过求解非线性方程组f(x)=0，得到系统的平衡点，并通过f(x)在平衡点的Jacobian矩阵特征值的情况判断平衡点的稳定性质。

2) 判断不稳定平衡点是否属于稳定平衡点边界上的点。

3) 利用稳定边界上的不稳定平衡点的稳定流形估计系统的稳定边界。

2.2 非线性映射的构造

设系统式(8)的平衡点满足双曲特性，那么存在同胚非线性映射[20]使非线性系统式(8)映射为线性系统：

(9)

式中：z与x具有相同的维数；Jr为非线性系统式(1)在不稳定平衡点的Jacobian矩阵A的Jordan标准型，即矩阵A与Jr满足：

(10)

式中：矩阵P为矩阵A的特征向量矩阵。

根据正规形理论，非线性系统式(8)映射为线性系统式(9)的具体过程如图1所示。

1) 将非线性系统式(8)在UEP进行Taylor级数展开，使其转变成为多项式形式的非线性系统:

(11)

2) 利用线性映射x=Py将系统式(11)映射为

(12)

3) 构造高阶多项式形式的非线性映射

(13)

通过非线性映射式(13)可将系统式(12)映射为线性系统式(9)。

图1 非线性映射流程图
Fig.1 Flow chart of nonlinear mapping

2.3 非线性映射的求解

下面利用待定系数法求解非线性映射式(13)的系数。

结合式(9)，取式(13)对时间t的导数，可得

(14)

将式(13)代入式(12)可得

(15)

联立式(14)和式(15)可得恒等式

(16)

令

(17)

(18)

(19)

根据式(19)可知式(13)的2次项系数a2可通过解线性方程组b2=0得到；得到2次项系数a2后系统的3次项系数a3可通过解线性方程组b3=0得到；以此类推便可得到式(13)的所有项的系数，其求解流程如图2所示。

图2 多项式系数求解流程图
Fig.2 Flow chart of solving polynomial coefficients

通过构造和求解非线性映射的过程可知，非线性映射具有统一的表达形式，且求解过程具有一定的程式化，便于计算机编程求解。

2.4 稳定边界的确定

考虑线性系统式(9)，由于状态矩阵为对角形式，因此每一特征值对应的特征向量Pi(i=1,2，…,N)有且仅有一个非零元素，即

Pi=[0 … 1(i)… 0]

(20)

设Jr具有负实部特征值为λ1,λ2，…,λN1，特征值对应的特征向量为ζ1,ζ2，…,ζN1；具有正实部特征值λ1,λ2，…,λN2，特征值对应的特征向量为η1,η2，…,ηN2，且N1+N2=N，那么线性系统式(2)的不变稳子空间由向量ζ1,ζ2,…,ζN1张成，不稳定子空间由η1,η2,…,ηN2张成[21]。根据不同特征值对应特征向量正交的特点可得特征向量ηj(j=1,2,…,N2)与稳定子空间的每一特征向量ζj(j=1,2,…,N1)正交。

(21)

因此，特征向量ηj与系统稳定子空间的向量z正交。

zTηj=0j=1,2,…,N2

(22)

结合式(20)可知式(22)等效为

zj=0j=1,2,…,N2

(23)

所以，线性系统的稳定子空间即稳定流形可表示为式(23)形式的N1维超曲面。因此结合式(13)、式(23)和线性映射x=Py可得非线性系统式(8)在不稳定平衡点的稳定流形为

(24)

根据2.1节流形理论估计系统稳定域的方法可知，通过对系统式(8)SEP稳定边界上所有UEP进行分析便可得到系统的稳定边界∂A

(25)

令w(x)满足：

(26)

则可通过边界参数w(x)w(xsep)判断系统状态是否为稳定平衡点吸引区内的状态。如果w(x)w(xsep)=0，则系统状态x位于稳定平衡点SEP稳定边界上；如果w(x)w(xsep)>0，则系统状态x位于稳定平衡点SEP吸引区内；如果w(x)w(xsep)<0，则系统状态x位于稳定平衡点SEP吸引区外。

3 着陆阶段结冰飞机的稳定域研究

3.1 高阶NF法估计稳定边界的有效性验证

高阶NF法以流形理论为基础，为说明高阶NF法的正确性，首先给出基于流形方法估计非线性系统稳定域的准确性，仿真结果如图3所示。

图3为基于流形方法得到的稳定边界与相同情况下系统真实边界的对比，用于证明流形方法估计系统稳定边界的有效性。其中彩色平面为流形法估计的稳定边界形状；图中红色点是通过仿真方法得到系统在特定范围内位于稳定域内或边界上的点，可准确反映出系统边界的位置和形状，可反映出系统稳定域的真实情况。图中所示稳定域内所用红点均位于流形法计算得到的稳定边界或包含于流形法的稳定边界内，因此可说明利用流形法估计稳定边界是有效的，且估计的边界具有较高精度。由于NF法估计的稳定域在多项式阶数达到一定水平时(有效阶数)将无限趋近于流形理论估计的稳定域，因此可以说明利用高阶NF法估计系统稳定域是有效的。下面以流形方法估计的稳定域作为参照，利用仿真的手段通过逐渐增加非线性映射阶数的方法确定利用NF法有效且能较为精确估计该统稳定域的阶数。仿真结果如图4和图5所示。

图3 流形方法有效性验证
Fig.3 Feasibility verification of manifold method

图4 NF法估计稳定域的有效阶数确定
Fig.4 Effective order of stability region estimation by NF method

图5 7阶NF法估计的稳定域
Fig.5 Stability region estimation by 7-order NF method

图4为NF法估计的稳定边界随多项式阶数变化情况，图5为7阶NF法估计的稳定边界与流形方法估计的稳定边界的对比。综合考虑图4和图5可知多项式阶数越高，NF法估计的稳定边界越接近真实情况，对于本文情况7阶NF法估计的稳定边界与系统真实边界情况基本吻合，因此后文将采用7阶NF法对结冰飞机的稳定边界变化及运动情况进行分析。

3.2 不同结冰程度对飞机稳定域的影响

首先计算干净飞机本体系统的平衡点并选取符合飞机飞行情况，且满足低速着陆状态要求的平衡状态作为工作点进行稳定域的分析，本文选择(47.00,0.156 5,0.170 5)为工作点，该点的特征值为

λ1,2=-0.599 7±2.249 0i

λ3=-0.097 8

通过计算发现此状态下系统的阻尼比ζ=0.257 7 较小，不满足国军标GJB185-86起降设计要求(该背景飞机采用放宽静稳定度技术)，因此利用式(7)对系统进行增稳控制补偿设计。通过设计可得控制参数为kα=-0.8，kq=-0.3，此时工作点特征值为

λ1,2=-1.500±1.823 7i

λ3=-0.039 6

此时阻尼比ζ=0.635 7满足设计标准。

其次，计算受控系统的平衡点并利用系统的Jacobian矩阵A的特征值判断平衡点的稳定性质。干净飞机在该工作点下存在两个平衡点A(47.00, 0.156 5, 0.170 5)和B(48.68, 0.125 2, 0.190 9)，通过上述方法可判断A为稳定平衡点SEP，B为不稳定平衡点UEP1，下面通过仿真的方法验证该UEP1为SEP边界上的不稳定平衡点，具体方法为对UEP1施加小扰动量Δα=0.05，如果系统运动稳定在SEP，则可证明UEP1为SEP边界上的不稳定平衡点，否则UEP1不在SEP边界上，仿真结果如图6所示。

图6为在UEP1添加扰动后系统状态x的动态响应图，可知添加扰动后系统动态x趋近于SEP，因此可以证明UEP1为SEP边界上的不稳定平衡点。利用同样的方法计算受控飞机在相同的工作点SEP，不同结冰程度下的平衡点，并判断平衡点的稳定性。通过计算发现系统的稳定平衡点SEP(47.00,0.156 5,0.170 5)保持不变，且在其边界上都仅存在一个不稳定平衡点，UEP计算结果如表2所示。

最后利用7阶NF方法求解结冰飞机在不同结冰情况的稳定边界，即UEP的稳定流形。根据本文2.3节和2.4节内容可得结冰飞机不同结冰程度的稳定边界表达式ws(xuep)，由于此过程通过计算机编程实现，故在此不进行详细说明。具体结果如图7和图8所示。

图6 扰动后的状态曲线
Fig.6 Dynamic response curves after disturbing

ηEPVt/(m·s-1)α/radq/(rad·s-1)0UEP148．680．12520．19090．07UEP248．400．12970．18750．12UEP348．200．13310．18510．20UEP447．880．13870．1814

图7 结冰因子对飞机稳定域的影响
Fig.7 Effect of icing factor on aircraft stability region

图8 结冰因子对UEP的影响
Fig.8 Effect of icing factor on UEP

图7为结冰飞机稳定域随结冰程度的变化情况，随着结冰程度的加剧，稳定域逐渐收缩，稳定边界下边缘逐渐接近平衡点，且结冰飞机稳定边界具有相似结构。图8为系统UEP随结冰程度的运动情况，随着结冰程度的加剧，UEP逐渐接近SEP，且结冰程度η=0.217 时，UEP与SEP重合，系统处于临界稳定状态，此时微小的扰动将会导致飞机的失稳；继续增加结冰程度(η>0.217)，此时系统SEP转变为UEP即系统稳定性质发生变化，而系统的稳定工作点沿速度减小方向逐渐远离原工作点，即系统工作状态发生漂移，此时飞机仍然处于临界稳定状态，但在此结冰程度情况下，飞机的飞行因丧失对原平衡状态的抗干扰能力而显得十分危险(将在后文进行详细说明)。

3.3 基于稳定域的结冰飞机动态性能分析

随着飞机结冰程度的加剧飞机的飞行稳定性将发生变化，因此根据结冰飞机飞行稳定性是否发生变化定义积冰对飞机影响的严重程度具有一定的意义，也从本质上反映出结冰的危害程度。

对于本文研究的飞行器而言定义飞机结冰后工作状态不发生漂移(SEP不变、飞行稳定性不变)的情况为轻度结冰情况，此时结冰因子η<0.217；定义结冰飞机工作状态发生漂移(SEP转变为UEP、飞行稳定性质发生变化)的情况为严重结冰情况，此时结冰因子η>0.217。下面利用高阶NF法给出两种结冰情况对飞行安全的威胁。

3.3.1 飞机轻度结冰情况

以着陆阶段飞机轻度结冰为例，分析轻度结冰的危害，选择结冰因子η=0.1进行分析。根据3.2节分析可知，飞机轻度结冰后，工作点SEP不会发生漂移，即飞机仍能保持原工作状态不变，但如果此时飞行员增加杆力使飞机快速改平，则飞机将面临失速危险，其飞行仿真结果如图9所示。

图9为着陆拉平减速阶段飞机轻度结冰后的运动情况。图中黄线为干净飞机的运动情况，而对于轻度结冰飞机而言，相同的操纵情况下使结冰飞机的运动情况首先沿黄线运动后继续沿红线运动，最终导致飞机失稳。原因是相同的操纵情况对于干净飞机，运动状态始终位于稳定域范围内，而对于轻度结冰(η=0.1)飞机，运动状态穿越其稳定域，飞机进入不稳定飞行状态，如不采取措施最终将导致飞行事故。综上所述，轻度结冰导致飞机非线性稳定域缩小，系统动态稳定性变差、抗干扰能力减弱、可操纵性变差。

下面给出高阶NF法中相关参数w(x)w(xsep)对飞行状态稳定性判断的方法。

图10为飞行过程中参数w(x)w(xsep)的变化情况，根据2.4节所述，对于干净飞机飞行过程中始终满足w(x)w(xsep)>0，因此干净飞机的飞行是安全稳定的；对于轻度结冰(η=0.1)飞机在临界稳定时间tc(critical time)以前参数w(x)w(xsep)>0，而当飞行时间t>tc时，参数w(x)w(xsep)<0，飞机进入不稳定飞行状态，如不采取措施最终将导致飞行事故，综上可知利用参数w(x)w(xsep)的飞行稳定性分析结果与利用NF法绘制稳定域的直观分析方法相同，而参数w(x)w(xsep)给出了运动状态与稳定边界之间的定量关系，分析过程不需借助图像，方法简单。

图9 轻度结冰后飞行动态
Fig.9 Flight dynamic of aircraft with mild icing

图10 基于NF法的运动稳定性分析
Fig.10 Stability analysis of movement by NF method

3.3.2 飞机严重结冰情况

以着陆阶段飞机严重结冰为例，分析严重结冰的危害，选择结冰程度η=0.25进行分析。根据3.2节分析可知，飞机严重结冰后，工作点SEP会发生漂移，原工作点处于临界稳定状态，即飞机很难维持原工作状态，飞机受微小扰动将进入失稳或趋近于另一平衡状态，但如果飞机经扰动后系统进入稳定工作状态，由于工作发生漂移，为保持原有工作状态，飞行员将进行相应的操纵，此时的操纵将成为飞机失稳的致命因素，上述情况的飞行仿真结果如图11所示。

图中UEP为原工作点，SEP为严重结冰后飞机的新稳定平衡点。蓝线为扰动后状态位于稳定域外时飞机的运动情况，由于状态位于稳定域之外，飞机在此后的运动将迅速失稳；绿线为扰动后状态位于稳定域内飞行员不进行操纵时飞机的运动情况，此时飞机将改变工作状态进入新的稳定情况；红线为扰动后状态位于稳定域内飞行员进行操纵时飞机的运动情况，由于新的状态变化较大，飞行员进行操纵以保持原飞行状态，此时飞机运动状态穿越稳定域，飞机将面临失速的危险。利用参数w(x)w(xsep)对严重结冰危害的分析过程与轻度结冰类似且同样可以得到上述结论，在此不进行重复说明。

图11 严重结冰情况后飞行动态
Fig.11 Flight dynamic of aircraft with severe icing

综上所述，在严重结冰情况下，飞机受到小扰动后，飞机将面临失速的危险，原因是严重结冰导致飞机在该工作点的稳定性质发生变化，恶化了系统的动态稳定性能，使系统丧失了在原平衡状态的抗干扰能力、且可操纵范围大大缩小。

4 结 论

本文基于高阶NF法研究了飞机的稳定域变化情况，并利用参数w(x)w(xsep)分析了不同结冰程度对飞机飞行安全产生的威胁，得到如下结论：

1) 高阶NF法可得到稳定边界的解析表达式，所估计的稳定边界具有较高精度，且适用于飞机结冰情况的稳定域研究。

2) 高阶NF法中参数w(x)w(xsep)可作为评估系统飞行安全性的指标，用于分析结冰程度对飞行运动的影响具有一定的优越性。

3) 不同结冰程度对飞机飞行稳定性质产生不同的影响，以此为依据定义轻度结冰和严重结冰情况，该定义方法从本质上反映出结冰对飞行安全的影响程度，因此该定义方法具有一定的科学性和可行性。轻度结冰导致飞机稳定域收缩、动态性能下降，抗干扰能力减弱；严重结冰情况导致飞机在平衡状态的飞行稳定性质发生变化，动态稳定性能恶化，丧失了对该工作点的抗干扰能力。

4) 飞机结冰对飞行安全产生严重威胁，尤其是在未察觉飞机结冰情况下，飞行员的操纵是加快结冰飞机失稳的主要原因和决定性因素。

[1] 屈亮， 李颖晖， 袁国强， 等. 基于相平面法的结冰飞机纵向非线性稳定域分析[J]. 航空学报， 2016, 37(3): 865-872.

QU L, LI Y H, YUAN G Q, et al. Longitudinal nonlinear stabilizing region for icing aircraft based on phase-plane method[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(3): 865-872 (in Chinese).

[2] HILTNER D W. A nonlinear aircraft simulation of ice contaminated tailplane stall[D]. Columbus: Ohio State University, 1998.

[3] Safety Advisor. Aircraft icing[EB/OL]. (2013-05-01)[2015-05-08]. http://www.aopa.org_media/-Files/AOPA/Home/Pilot%20Resources/ASI/Safety%20Advisors/sall.pdf.

[4] HAROLD E. ADDY J R. Ice Accretions and icing effects for modern airfoils: NASA/TP-2000-210031[R]. Wa-shington, D.C: NASA, 2000.

[5] MILLER R, RIBEENS W. The Effects of icing on the longitudinal dynamics of an icing research aircraft: AIAA-1999-0636[R]. Reston: AIAA, 1999.

[6] BRAGG M B, HUTCHISON T, MERRET J, et al. Effect of ice accretion on aircraft flight dynamic: AIAA-2000-0360[R]. Reston: AIAA, 2000.

[7] 裘燮纲, 韩凤华. 飞机防冰系统[M]. 北京: 航空专业教程编审组, 1985.

QIU X G, HAN F H. Aircraft anti-icing system[M]. Beijing: Aviation Professional Editors Publishing Group, 1985 (in Chinese).

[8] 王明丰, 王立新, 黄成涛. 积冰对飞机纵向操稳特性的量化影响[J]. 北京航空航天大学学报, 2008, 34(5): 592-595.

WANG M F, WANG L X, HUANG C T. Computational effects of ice accretion on aircraft longitudinal stability and control[J]. Journal of Beijing University of Aeronautics and Astronautics, 2008, 34(5): 592-595 (in Chinese).

[9] 徐忠达, 苏媛, 曹义华. 积冰对飞机操纵性的影响与仿真[J]. 北京航空航天大学学报, 2012, 38(7): 941-946.

XU Z D, SU Y, CAO Y H. Simulation of ice effects on aircraft controllability[J]. Journal of Beijing University of Aeronautics and Astronautics, 2012, 38(7): 941-946 (in Chinese).

[10] 周莉, 徐浩军, 杨哲. 飞机在结冰条件下的最优边界保护方法[J]. 上海交通大学学报, 2013, 47(8): 1217-1221.

ZHOU L, XU H J, YANG Z. Optimal boundary protection method for aircraft under icing conditions[J]. Journal of Shanghai Jiao Tong University, 2013, 47(8): 1217-1221 (in Chinese).

[11] POKHARIYAL D, BRAGG M B, HUTCHISON T, et al. Aircraft flight dynamics with simulated ice accretion: AIAA-2001-0541[R]. Reston: AIAA, 2001.

[12] 高浩, 周志强. 高机动飞机迎角全局稳定性研究[J]. 航空学报, 1987, 8(11): 562-571.

GAO H, ZHOU Z Q. A study of the global stability of high performance aircrafts at high angle-of-attack[J]. Acta Aeronautica et Astronautica Sinica, 1987, 8(11): 562-571 (in Chinese).

[13] 何植岱, 郭文. 非线性飞行稳定性研究的新综合方法[J]. 空气动力学学报, 1990, 8(2): 143-151.

HE Z D, GUO W. A novel synthetical method for study nonlinear flight stability[J]. Acta Aerodynamics Sinica, 1990, 8(2): 143-151 (in Chinese).

[14] ROBERT C, ALLEN, HARRY G K. Safe set protection and restoration for unimpaired and impaired aircraft: AIAA-2012-4822[R]. Reston: AIAA, 2012.

[15] 王爽， 詹浩. 飞行最大可控边界集及其机动边界保护控制[J]. 西北工业大学学报, 2014, 32(4): 523-528.

WANG S, ZHAN H. The safe-set of aircraft and maneuverability envelope protection[J]. Journal of Northwestern Polytechnical University, 2014, 32(4): 523-528 (in Chinese).

[16] 李颖晖, 张保会, 李勐. 电力系统稳定边界的研究[J]. 中国电机工程学报, 2002, 22(3): 72-77.

LI Y H, ZHANG B H, LI M. Study on electrical power system stability boundary[J]. Proceedings of the CSEE, 2002, 22(3): 72-77 (in Chinese).

[17] SAHA S, FOUAD A A, KLIEMANN W H, et al. Stability boundary approximation of a power system using the real normal form of vector fields[J]. IEEE Transactions on Power Systems, 1997, 12(2): 797-802.

[18] CHIANGH D, HIRSCH M W, WU F F. Stability region of nonlinear autonomous dynamical system[J]. IEEE Transactions on Automatic Control, 1988, 33(1): 16-27.

[19] CHIANG H D, JAMES S, THORP. Stability regions of nonlinear dynamical system: A constructive methodology[J]. IEEE Transactions on Automatic Control, 1988, 34(12): 1229-1241.

[20] WANG D. Anintroduction to the normal form theory of ordinary differential equation[J]. Advances in Mathematic, 1990, 86(1): 38-71.

[21] 李颖晖, 张保会. 正规形理论在电力系统稳定性研究中的应用(二)——电力系统主导不稳定平衡点上局部流形的计算[J]. 电力自动化设备, 2003, 23(7): 1-4.

LI Y H, ZHANG B H. Application of normal form in study of power system stability part 2: Calculation of local manifolds on controlling unstable equilibrium point of electric power system[J]. Electric Power Automation Equipment, 2003, 23(7): 1-4 (in Chinese).

(责任编辑： 李明敏)

URL：www.cnki.net/kcms/detail/11.1929.V.20161103.1633.004.html

Nonlinearstabilityregionoficingaircraftduringlandingphasebasedonnormalformmethod

ZHENGWuji,LIYinghui*,QULiang,XUHaojun,YUANGuoqiang

AeronauticsandAstronauticsEngineeringCollege,AirForceEngineeringUniversity,Xi’an710038,China

Icingcancauseflightenvelopeshrinkandthusposesagreatthreattoflightsafety;therefore,researchonnonlinearstabilityregionofanaircraftaftericingissignificantlyimportantforthedesignofflightenvelopeprotectionsystemandtheimprovementflightsafety.Atransportistakenasanobjectofthisstudy.Takingthenonlinearaerodynamiccharacteristicsintoaccount,thelongitudinalnonlineardynamicmodelwithstabilityaugmentationcontrolisobtained.Basedonthetheoryofmanifoldandnormalform,thelongitudinalnonlinearstabilityboundaryanditsanalyticexpressionareobtained.Viadynamicsimulation,thestabilityboundarybasedonnormalformisjustifiedtobefeasibleandaccurate.Finally,icingfactorimpactingonaircraftstabilityregionandthemechanismforaccidentsofanicingaircraftisstudiedatthelandingphase.Resultsshowthatmildicingconditioncancausetheshrinkofnonlinearstabilityregion,whilesevereicingconditioncanchangetheaircraftstability.Whenicingofanaircraftisnotdetected,regularmanipulationcanresultinflightaccidentastheaircraftstatecanbeoutoftheicingnonlinearstabilityregion.Theresultscanprovidesomereferenceforflightenvelopeprotectionundericingcondition.

aircrafticing;stabilityregion;nonlinear;normalformtheory;transportaircraft

2016-08-25；Revised2016-09-20；Accepted2016-10-25；Publishedonline2016-11-031633

s：NationalNaturalScienceFoundationofChina(61374145);NationalBasicResearchProgramofChina(2015CB755805)

.E-mailliyinghui66@163.com

2016-08-25；退修日期2016-09-20；录用日期2016-10-25； < class="emphasis_bold">网络出版时间

时间：2016-11-031633

www.cnki.net/kcms/detail/11.1929.V.20161103.1633.004.html

国家自然科学基金 (61374145)； 国家“973”计划 (2015CB755805)

.E-mailliyinghui66@163.com

郑无计， 李颖晖， 屈亮， 等． 基于正规形法的结冰飞机着陆阶段非线性稳定域J． 航空学报,2017,38(2):520714.ZHENGWJ,LIYH,QUL,etal.NonlinearstabilityregionoficingaircraftduringlandingphasebasedonnormalformmethodJ.ActaAeronauticaetAstronauticaSinica,2017,38(2):520714.

http://hkxb.buaa.edu.cnhkxb@buaa.edu.cn

10.7527/S1000-6893.2016.0279

V328

A

1000-6893(2017)02-520714-11