对先进战斗机而言，在大迎角范围内能够进行大幅度可控的快速机动是一项必备的能力。当战斗机进行过失速机动时，绕飞机的流场会出现分离涡等复杂流动现象，此时作用在飞机上的气动力(矩)将呈现明显的非线性和非定常特性，即气动力(矩)不仅取决于当前时刻的运动状态，还与运动的时间历程有关[1]。对于这种情况，传统的气动导数模型已经很难适用[2]，需要探索更加准确高效的非线性非定常气动力建模方法。
相关建模方法已经有不少研究[3-5]，常见的包括阶跃响应模型[1]、Volterra级数模型[6-7]、状态空间模型[8-9]、微分方程模型[10-12]等。随着机器学习的蓬勃发展，近年来出现了将机器学习应用到气动力建模的尝试，如神经网络[13]、模糊逻辑[14]、支持向量机[15-16]等，这类方法的优点在于将系统视为黑箱，利用机器学习对复杂系统的描述能力预测输出。其中，支持向量机(Support Vector Machine, SVM)是Vapnik[17]在1995年基于统计学习理论提出的一种针对小样本情况的机器学习方法，其最初应用于模式识别，后来逐渐扩展到系统辨识、智能控制等领域[18-19]。最小二乘支持向量机(Least Squares Support Vector Machine, LS-SVM)[20]是标准SVM的一种改进方法，其将二次规划问题转化为线性方程组求解，降低了计算复杂度。
在现有的非定常气动力建模过程中，大多是以模型在风洞中的强迫振动试验所得到的数据为初始样本进行，其流程如图 1所示。在某一减缩频率下，当试验对象经过多个周期的强迫振动之后，试验对象和其周围的流动以及产生的气动力均呈现出稳定的周期性变化特征，此时的运动是一个稳定的动态过程，此过程的气动力数据即是该减缩频率下的稳定滞环数据，采集该数据作为样本。接下来，在另一减缩频率下继续进行试验。当所有试验完成之后，获得多个减缩频率下的稳定滞环数据。最后，基于这些数据进行建模。建模完成之后，模型被用来预测其他减缩频率下的气动力响应。但是，这些方法针对的对象始终是稳定的滞环，对于气动力是如何达到稳定滞环的，也就是气动力的初始发展过程，即图 1中气动力稳定之前的阶段没有考虑。飞机的快速机动过程是一个动态的变化过程，可以认为在这个过程中飞机不会持续运动到进入稳定滞环，而是处于初始的非稳定阶段。因此，初始非稳定阶段的气动力预测相比稳定段，更具有实际意义。这里以某机翼为例，由CFD计算的在不同减缩频率k下做俯仰运动时的升力系数CL滞环，如图 2所示，初始迎角α为30°。图中，实线表示稳定滞环，虚线表示进入稳定滞环的初始非稳定过程。可以发现，气动力明显存在一个进入稳定滞环的过程。所以，仅仅预测稳定滞环是不够的，更重要的是能够预测初始非稳定过程的气动力。
|
| 图 1 基于风洞试验数据的气动力建模方法示意图 Fig. 1 Schematic diagram of aerodynamic modeling method based on wind tunnel test data |
|
| 图 2 某机翼做俯仰运动时的升力系数滞环(Ma=0.4, Re=3.0×106) Fig. 2 Lift coefficient hysteresis of a wing for pitching motion at Ma=0.4 and Re=3.0×106 |
传统方法基于多个减缩频率的稳定滞环数据进行建模，存在的另一个问题是，系统每次仅受到单个频率的激励，只能激发出系统的局部线性特性，而大迎角范围内的气动力是非线性的。对于战斗机的过失速机动飞行而言，真实的飞行过程以纵向俯仰为例，战斗机在拉升的过程中快速进入过失速机动，该过程几乎不可能是标准的振动，更不可能用单一频率的稳定振动进行模拟。为了能够正确地模拟战斗机动态进入过失速机动过程以及产生的非定常气动力(矩)，不仅需要能够包含非稳定动态特性影响，同时需要反映不同频率的非线性耦合影响。
本文根据Morelli等[21-22]考虑多操纵面非线性耦合影响而提出的正交相位优化多正弦激励方法，并合理选择输入形式，设计了用于构建大迎角俯仰动态过程的非线性激励模型。以某机翼为例，通过CFD技术获取了激励模型下的气动响应，并以此为构建大迎角非定常气动力模型的样本数据，利用LS-SVM方法建立了仅通过一次激励，就能预测在幅值和频率范围内任意运动的非线性非定常气动力模型。算例结果表明，该模型不仅能够正确地预测不同频率下的稳定滞环，同时能预测初始的非稳定过程。
1 激励模型 1.1 非稳定动态过程和非线性系统辨识首先借助振动理论讨论一般线性系统在外部激励下的响应，分析非定常气动力响应的变化过程。以一个典型的二阶振动器系统为例，其运动方程为
| $ \ddot y + 2\zeta p\dot y + {p^2}y = {Y_0}{p^2}{\rm{cos}}(\omega t) $ | （1） |
式中：y为位移；ζ为阻尼比；p为自然频率；Y0p2cos(ωt)为外部输入。
基于振动理论，位移响应y包含两部分，分别为自由振动y1和强迫振动y2，即
| $ {y = {y_1} + {y_2}} $ | （2） |
| $ {{y_1} = {{\rm{e}}^{ - \zeta pt}}(A{\rm{sin}}(qt) + B{\rm{cos}}(qt))} $ | （3） |
| $ {{y_2} = Y{\rm{cos}}(\omega t - \varphi )} $ | （4） |
式中：A、B为系数，具体表达式可见振动理论相关教材；
对于一个稳定系统，自由振动y1会随时间的推移而逐渐衰减到零，最后仅剩下强迫振动y2。以某个振荡器[23]为例，不考虑非线性环节时，在外部输入为正弦的情况下，系统的自由振动响应和强迫振动响应如图 3所示。可见，系统响应在开始的时候，存在一个由自由振动响应和强迫振动响应共同构成的非稳定段，这个过程是一个非稳定动态过程。随着时间的推移，自由振动衰减到零，只剩下强迫振动，此时的响应达到稳定，这个过程是一个稳定动态过程。传统的非定常气动力建模方法，仅采用了稳定段的非定常气动力数据，而忽略了初始的非稳定过程。换言之，只有强迫振动数据被采集和使用，而自由振动数据被遗漏了，这将导致对系统响应变化过程的信息缺失。因此，基于这些数据的模型只能预测最终的稳定滞环，而不能预测进入滞环的初始非稳定过程。
其次，大迎角区的气动力会呈现明显的非线性特征。对于非线性系统而言，单个频率激励只能激发出系统在该频率下的局部特性，无法激发出系统中不同频率耦合作用的影响。因此，如果仅分别使用不同频率的稳定滞环为样本数据进行大迎角非定常气动力建模，所获得的模型只能够反映不同频率下稳定的非定常气动力的迟滞环特性，不能反映系统中频率耦合的非线性影响。
因此，为了建立能够正确计算过失速机动过程的非定常气动力模型，首先需要构建能够激发非线性动态过程中频率耦合作用的激励模型，其既包含了非稳定动态过程的影响，又适于非线性系统辨识。
1.2 激励的选择和设计对于系统辨识和机器学习问题，激励的选择和设计至关重要，直接影响了模型的预测精度。从系统辨识的角度，选择训练数据的关键是充分激发出系统特性。对于线性系统辨识，有多种输入，例如扫频[24]、3211[25]等可以采用。对于非线性系统辨识，缺乏通用的输入信号。
根据Prazenica和Kurdila[26]的研究，在频域表示线性系统特性的频率响应函数是一阶的，可以表示为H1(jω)，其中ω是频率，j是虚部。以扫频激励为例，其触发到了不同频率下的系统特性，因此适用于辨识一阶频率响应函数H1(jω)，即适用于辨识线性系统。表示非线性系统特性的频率响应函数是高阶的，包含了H1(jω)、H2(jω1, jω2)、H3(jω1, jω2, jω3)等。换言之，非线性系统的响应不仅包含了单一频率的影响，还包含了不同频率的耦合影响，此时扫频激励不再适用。因此，激励的选择和设计必须能够反映不同频率的耦合作用。
这里采用Morelli等[21-22]设计的用于飞机多操纵面耦合激励的输入——正交相位优化多正弦，作为训练输入。为了对飞机多个操纵面同时进行激励，以实现飞机多通道气动导数的同时辨识，Morelli等需要使得不同操纵面的激励在时域和频域上都互相正交，以实现互不干扰。这种输入通过在不同通道将多个互相正交的正弦信号叠加，然后利用优化方法调节相位以避免不同分量相加使得幅值过大，导致系统受到激励后偏离参考点过远，最终使得不同通道之间的激励互相正交，同一通道中不同分量之间也互相正交，并覆盖感兴趣的频段。借鉴这种思路，这里将多通道合并到一个通道，同时仍然保持正交性，构造了一种能体现不同频率耦合作用的输入，用于非线性系统辨识。本文作者[27]成功将这种输入应用到Volterra级数高阶核的辨识，实现了翼型在跨声速下小迎角范围内的非线性非定常气动力建模。
正交相位优化多正弦输入的表达式为
| $ \mathit{\boldsymbol{u}} = \sum\limits_{i = 1}^M {{A_i}} {\rm{sin}}\left( {\frac{{2\pi it}}{T} + {\varphi _i}} \right) $ | （5） |
式中：u为输入；Ai为幅值；T为输入时长；φi为相位；M为选择的正弦总数。
根据Schroeder[28]的研究，激励信号要尽可能充分的激励系统，需要激励信号的能量尽可能大，衡量信号能量的是均方根(Root Mean Square, RMS)，表示为
| $ {\rm{RMS}} (\mathit{\boldsymbol{u}}) = \sqrt {\frac{{{\mathit{\boldsymbol{u}}^{\rm{T}}}\mathit{\boldsymbol{u}}}}{N}} $ | （6） |
式中：N为输入的离散点数。同时，还要使得系统响应不偏离参考点过远，保证这一点需要激励的幅值差值尽可能小，即max(u)－min(u)尽可能小。综合考虑这两个因素，Schroeder提出，可以利用相对峰值因子(Relative Peak Factor, RPF)最小来确定输入相位，其定义为
| $ {\rm{RPF}} (\mathit{\boldsymbol{u}}) = \frac{{{\rm{max}}(\mathit{\boldsymbol{u}}) - {\rm{min}}(\mathit{\boldsymbol{u}})}}{{2\sqrt 2 \sqrt {({\mathit{\boldsymbol{u}}^{\rm{T}}}\mathit{\boldsymbol{u}})/N} }} $ | （7） |
该输入设计的具体过程如下：
步骤1 根据研究对象的实际情况确定输入的时长T和频带[fmin, fmax]，一般要求fmin≥2/T，由此可以得到频率精度Δf=1/T和M=fix[(fmax－fmin)/Δf]+1，其中fix(·)表示向下取整。
步骤2 根据需要确定功率谱密度，这里简单取均匀功率谱，选择式(5)中各个分量的幅值相等，即
步骤3 随机给定各个分量的相位初值φi∈(－π, π]，然后运用优化方法调节相位，使得各分量相加后的输入的相对峰值因子最小，最后得到需要的相位分量φi。
步骤4 根据前3步确定的参数构成的输入作小幅相位偏移，使得输入从零开始，到零终止，便于对系统施加激励，最终得到需要的输入。
本文算例所用的正交相位优化多正弦输入，取时长T=50 s，频率范围为fmin=2/50 Hz，fmax= 2 Hz，
|
| 图 4 正交相位优化多正弦输入 Fig. 4 Orthogonal phase-optimized multisine input |
在训练之前，必须根据实际问题的物理背景选择合理的输入形式。因为非定常气动力不仅依赖于当前的运动状态，还依赖于之前的运动过程，所以输入形式的选择必须反映这一点。在已有的利用机器学习进行气动力建模的研究中，减缩频率或者人为定义的等效减缩频率常常作为输入参数之一。例如，Ignatyev[29]、Wang[15]和史志伟[30]等，所用训练数据均来自于不同减缩频率下的强迫振动试验。然而，如前所述，真实的飞行过程几乎不可能是标准的振动，因而不存在固定频率。建立一个不依赖于频率且具有通用性的气动力模型是本文研究的问题之一。
这里选择的输入形式为
| $ \mathit{\boldsymbol{x}} = {[\alpha (\tau ),\alpha (\tau - 1), \cdots ,\alpha (\tau - m)]^{\rm{T}}} $ | （8） |
式中：α为迎角；m为当前时刻的气动力对之前运动过程的依赖长度。此处没有引入减缩频率，因为在激励信号的设计中已经通过频带的选择考虑了频率的影响。
2 最小二乘支持向量机 2.1 基本理论给定训练样本集(xi, yi)，i=1, 2, …, n，xi∈ Rd，yi∈ R，SVM回归的基本思想是将样本所在的输入空间Rd通过非线性映射φ映射到特征空间φ(x)，然后在该空间构造回归函数y=f(x)wTφ(x)+b，其中w为加权系数，b为偏差。参数(w, b)通过结构风险最小化原则确定，即求解如下优化问题：
| $ \left\{ {\begin{array}{*{20}{l}} {\mathop {{\rm{min}}}\limits_{\mathit{\boldsymbol{w}},b,\xi } }&{R = \frac{1}{2}{{\left\| \mathit{\boldsymbol{w}} \right\|}^2} + \frac{1}{2}c\sum\limits_{i = 1}^n {\xi _i^2} }\\ {{\rm{s}}{\rm{.t}}{\rm{.}}}&{{y_i} = {\mathit{\boldsymbol{w}}^{\rm{T}}}\mathit{\boldsymbol{\varphi }}({\mathit{\boldsymbol{x}}_i}) + b + {\xi _i}\;\:i = 1,2, \cdots ,n} \end{array}} \right. $ | （9） |
式中：R为结构风险损失函数；ξi为预测误差；
传统的学习方法主要基于经验风险最小化原则，即寻求样本集上的经验风险最小，其理论基础是当训练用的样本数据趋于无穷时，经验风险收敛于真实风险。但在实际问题中样本数一般是有限的，此时经验风险最小并不能保证真实风险也达到最小，因此基于该原则的诸多算法如神经网络常出现“过学习”问题。根据统计学习理论，学习机器的真实风险由经验风险和置信范围两部分构成，要降低真实风险，除了要降低经验风险，还要缩小置信范围。置信范围与学习机器的复杂度和训练样本数有关，在样本数不变的情况下，缩小置信范围就需要降低复杂度。减小|| w||能降低复杂度，因此将||w||2/2作为式(9)中的一部分。
用Lagrangian方法求解式(9)的等式约束优化问题，得
| $ \begin{array}{*{20}{c}} {L(\mathit{\boldsymbol{w}},b,{\xi _i},{\lambda _i}) = \frac{1}{2}{{\left\| \mathit{\boldsymbol{w}} \right\|}^2} + \frac{1}{2}c\sum\limits_{i = 1}^n {\xi _i^2} - }\\ {\sum\limits_{i = 1}^n {{\lambda _i}} ({\mathit{\boldsymbol{w}}^{\rm{T}}}\mathit{\boldsymbol{\varphi }}({\mathit{\boldsymbol{x}}_i}) + b + {\xi _i} - {y_i})} \end{array} $ | （10） |
式中：λi为Lagrangian乘子。最优解为
| $ \frac{{\partial L}}{{\partial \mathit{\boldsymbol{w}}}} = 0,\frac{{\partial L}}{{\partial b}} = 0,\frac{{\partial L}}{{\partial {\xi _i}}} = 0,\frac{{\partial L}}{{\partial {\lambda _i}}} = 0 $ | （11） |
对应条件为
| $ {\mathit{\boldsymbol{w}} = \sum\limits_{i = 1}^n {{\lambda _i}} \mathit{\boldsymbol{\varphi }}({\mathit{\boldsymbol{x}}_i})} $ | （12） |
| $ {\sum\limits_{i = 1}^n {{\lambda _i}} = 0} $ | （13） |
| $ {c{\xi _i} = {\lambda _i}} $ | （14） |
| $ {{y_i} = {\mathit{\boldsymbol{w}}^{\rm{T}}}\mathit{\boldsymbol{\varphi }}({\mathit{\boldsymbol{x}}_i}) + b + {\xi _i}} $ | （15） |
将式(12)和式(14)代入式(15)，消去w和ξi，得
| $ {y_i} = \sum\limits_{l = 1}^n {{\lambda _l}} \langle \mathit{\boldsymbol{\varphi }}({\mathit{\boldsymbol{x}}_i}),\mathit{\boldsymbol{\varphi }}({\mathit{\boldsymbol{x}}_l})\rangle + b + \frac{1}{c}{\lambda _i} $ | （16） |
式中：〈· 〉为内积。定义核函数K(xi, xl)=〈φ(xi), φ(xl)〉，将式(13)和式(16)合并，得方程组的矩阵形式为
| $ \left[ {\begin{array}{*{20}{c}} {\mathit{\boldsymbol{K}} + {c^{ - 1}}{{\bf{1}}_n}}&{{{\bf{1}}_n}}\\ {{\bf{1}}_n^{\rm{T}}}&{\bf{0}} \end{array}} \right]\left[ {\begin{array}{*{20}{l}} \lambda \\ b \end{array}} \right] = \left[ {\begin{array}{*{20}{l}} \mathit{\boldsymbol{y}}\\ 0 \end{array}} \right] $ | （17） |
式中：Kil=K(xi, xl)；λ=[λ1, λ2, …, λn]T；y=[y1, y2, …, yn]T；1n为n个1的列向量。
为了求解式(17)以得到参数(λ, b)，常用方法是利用核函数矩阵K的对称正定性质对方程组做变换，然后利用Cholesky分解求解，具体方法可见文献[31]。已知(λ, b)，即得到了所需的LS-SVM模型：
| $ f(\mathit{\boldsymbol{x}}) = \sum\limits_{i = 1}^n {{\lambda _i}} K(\mathit{\boldsymbol{x}},{\mathit{\boldsymbol{x}}_i}) + b $ | （18） |
SVM通过引入核函数避免了在高维特征空间直接计算内积，极大简化了计算。根据泛函的有关理论，任意满足Mercer条件的对称函数都能作为核函数，最常采用的是高斯核函数：
| $ K({\mathit{\boldsymbol{x}}_i},{\mathit{\boldsymbol{x}}_l}) = {\rm{exp}}\left( { - \frac{{{{\left\| {{\mathit{\boldsymbol{x}}_i},{\mathit{\boldsymbol{x}}_l}} \right\|}^2}}}{{2{\sigma ^2}}}} \right) $ | （19） |
核函数选定后，模型参数的确定主要是核参数σ和正则化参数c的确定。某组参数(σ, c)所对应的泛化误差常采用留一交叉验证法(Leave-One-Out Cross-Validation)计算，即给定一组(σ, c)之后，每次将训练样本中的一个样本作为测试集，剩下的作为训练集，计算预测误差，如此重复，直到每个样本都作为过测试集，最后将所有预测误差的平均值作为该组参数下的泛化误差。如何找到使得泛化误差最小的参数(σ, c)，可采用MATLAB优化工具箱中的fmincon函数，其利用共轭梯度方法进行搜索[31]。
2.3 训练方法输入信号选择以后，接下来利用CFD计算气动力。然后根据选择的输入形式构造输入向量xi，输出分别为升力系数、阻力系数和俯仰力矩系数。最后对LS-SVM进行训练。具体步骤如下：
步骤1 按照1.2节构造输入，然后由CFD计算气动力响应。
步骤2 根据1.3节构造输入向量xi，输出yi为气动力响应，i=1, 2, …, n。
步骤3 给定模型参数的初值，包括式(9)中的c和式(19)中的σ。
步骤4 利用MATLAB的fmincon函数，以留一交叉验证法计算的平均预测误差为优化目标，寻找最优的模型参数(σ, c)。
步骤5 回到步骤3和步骤4，可重复多次，避免(σ, c)的局部最优，最终选择使得平均预测误差最小的模型参数(σ, c)。
步骤6 求解式(17)，得到参数(λ, b)，即完成了模型的训练。在已知测试数据后，由式(18)即可计算模型的预测结果。
3 训练以某机翼在亚声速下大迎角范围内做俯仰运动为例，验证LS-SVM模型预测非线性非定常气动力的能力和所采用设计输入的有效性。训练所用数据由CFD计算得到，与所用模型无关。计算条件为Ma=0.4，Re=3.0×106，采用刚性动网格。机翼外形参数如图 5所示，剖面翼型为NACA-64A204，网格图如图 6所示。
|
| 图 5 机翼外形图 Fig. 5 Wing shape |
|
| 图 6 机翼网格生成图 Fig. 6 Sample grid generated around wing |
本文所用CFD求解器是NASA的Langley研究中心开发的开源计算工具CFL3D[32]，其采用有限体积法和Spalart-Allmaras湍流模型求解雷诺平均Navier-Stokes (RANS)方程。CFL3D已经在航空领域通过多个算例得到了广泛验证[33-35]。由CFL3D计算的机翼在当前状态下的定常气动力如图 7所示。
|
| 图 7 机翼在Ma=0.4、Re=3.0×106下的定常气动力系数 Fig. 7 Steady aerodynamic coefficients of wing at Ma=0.4, Re=3.0×106 |
以设计的输入由CFD进行非定常计算，得到气动力数据，然后进行训练。训练完成后，以训练数据作为测试样本检验训练效果，图 8和图 9分别为以不同基准迎角做俯仰运动下，升力系数、阻力系数和俯仰力矩系数的LS-SVM预测结果和CFD计算结果。Case 1的基准迎角为0°，幅值为80°；Case 2的基准迎角为30°，幅值为30°。关于式(8)中时间依赖长度m的选择，需要尝试不同的时间长度，然后根据训练误差进行选择，这里Case 1取m=10，Case 2取m=30。由图 8和图 9可见预测结果和CFD结果非常吻合，验证了LS-SVM对训练数据的学习能力。
|
| 图 8 Case 1：初始迎角α0=0°时LS-SVM预测结果与训练数据对比(前10 s) Fig. 8 Case 1: LS-SVM predictions and training data when initial angle of attack α0=0°(first 10 s) |
|
| 图 9 Case 2：初始迎角α0=30°时LS-SVM预测结果与训练数据对比(前10 s) Fig. 9 Case 2: LS-SVM predictions and training data when initial angle of attack α0=30° (first 10 s) |
模型训练完成后，需要验证其在其他输入下的预测准确性，也即验证其推广能力。这里分别对机翼在不同基准状态下俯仰做正弦和扫频运动时的气动力进行计算验证。
4.1 Case 1Case 1的基准迎角为0°，幅值为80°，机翼做俯仰运动时迎角变化为
| $ \left\{ {\begin{array}{*{20}{l}} {\alpha (t) = 80{\rm{sin}}(1.191t)}\\ {\alpha (t) = 80{\rm{sin}}(3.573t)}\\ {\alpha (t) = 80{\rm{sin}}(7.147t)}\\ {\alpha (t) = 80{\rm{sin}}(0.142{t^2})} \end{array}} \right. $ | （20） |
前3种正弦运动对应的减缩频率k分别为0.01、0.03、0.06，第4种为扫频运动，对应的减缩频率从0线性增加到0.08。
由LS-SVM计算的测试运动的气动力响应和CFD验证结果对比如图 10所示。3种正弦运动的气动力滞环包括初始非稳定动态过程，都和CFD结果一致。图中虚线加粗部分是进入滞环的非稳定动态过程，可以发现，当减缩频率为0.01和0.03时，气动力的非稳定段与稳定段几乎重合，差异不明显；当减缩频率增加到0.06时，气动力的非稳定段与稳定段出现了差异。通过对比可以看到，不仅稳定滞环部分得到了准确预测，进入稳定滞环的初始非稳定过程也得到了准确预测。
|
| 图 10 初始迎角α0=0°时LS-SVM预测结果与CFD结果 Fig. 10 Results of LS-SVM predictions and CFD data when initial angle of attack α0=0° |
除了3种正弦运动，扫频运动的预测结果与CFD结果也吻合很好，表明通过采用文中的训练数据设计方法和输入形式，避免了以减缩频率为参数，不仅能预测正弦运动，理论上能预测范围内的任意运动，使得模型更具有普适性。
4.2 Case 2Case 2的基准迎角为30°，幅值为30°，机翼做俯仰运动时迎角变化为
| $ \left\{ {\begin{array}{*{20}{l}} {\alpha (t) = 30 + 30{\rm{sin}}(1.191t)}\\ {\alpha (t) = 30 + 30{\rm{sin}}(3.573t)}\\ {\alpha (t) = 30 + 30{\rm{sin}}(7.147t)}\\ {\alpha (t) = 30 + 30{\rm{sin}}(0.142{t^2})} \end{array}} \right. $ | （21） |
前3种正弦运动的减缩频率分别为0.01、0.03和0.06，第4种为扫频运动，对应的减缩频率从0线性增加到0.08。
由LS-SVM计算的测试运动的气动力响应和CFD验证结果对比如图 11所示。3种正弦运动中的气动力滞环随着减缩频率的增大，非稳定动态过程表现逐渐明显。包括非稳定动态过程在内的整个滞环，都实现了准确预测，即使升力系数、阻力系数和俯仰力矩系数随减缩频率的变化呈现不同的特征，预测结果与CFD结果都非常吻合。
|
| 图 11 初始迎角α0=30°时LS-SVM预测结果与CFD结果 Fig. 11 Results of LS-SVM predictions and CFD data when initial angle of attack α0=30° |
与Case 1类似，扫频运动的预测结果与CFD结果也吻合很好，再次表明通过采用文中的训练数据设计方法和输入的形式选择，不仅能预测正弦运动，理论上能预测幅值和频率范围内的任意运动，使得模型更具有普适性。
Case 2相比Case 1，基准迎角发生了变化，气动力的初始非稳定过程表现得更为明显，也就是说，在Case 2状态下，气动力的初始段和稳定滞环存在明显差异。根据前面的分析，飞机的真实机动过程始终处于初始段，不会到达稳定滞环。因此可以认为，在当前状态下，仅依靠稳定滞环来预测飞机机动过程中的非定常气动力将存在很大误差。
4.3 以强迫振动为训练输入时的预测现在将1.1节讨论的强迫振动作为训练输入，其他不变，选用式(21)中的第2个正弦运动为算例，验证气动力的预测效果。由CFD计算的升力系数的前3个周期响应如图 12所示。可以看到，气动力很快趋于周期性变化。因此，可以认为后两个周期的气动力响应是稳定的。分别将后两个周期的稳定振动数据和包含初始非稳定段在内的全部3个周期的振动数据作为训练数据，然后预测同样正弦运动下的升力系数响应。预测结果的时间响应曲线(前0.6 s)和滞环曲线如图 13所示。标注“稳定输入”的曲线为采用后两个周期的稳定数据为训练数据所得到的预测结果，标注“全部输入”的曲线为采用全部3个周期的数据为训练数据所得到的预测结果。
|
| 图 12 初始迎角α0=30°时由CFD计算的升力系数训练输入(减缩频率k=0.03) Fig. 12 Training data of lift coefficient from CFD with a reduced frequency of k=0.03 when initial angle of attack α0=30° |
|
| 图 13 初始迎角α0=30°时不同训练输入下LS-SVM预测结果与CFD结果(减缩频率k=0.03) Fig. 13 Results of LS-SVM predictions and CFD data with a reduced frequency of k=0.03 and different training inputs when initial angle of attack α0=30° |
通过和CFD结果对比可以发现，采用稳定振动数据进行训练，可以准确预测稳定滞环，但进入稳定滞环的初始非稳定过程无法实现预测。采用包含初始变化过程的振动数据进行训练之后，可以实现气动力整个发展过程的准确预测。这个例子再次表明，采用强迫振动输入只能预测到稳定滞环，只有当包含初始的自由振动之后，进入滞环的非稳定过程才能得到预测。因此，基于强迫振动试验数据的建模方法难以用于飞机快速机动过程中非定常气动力的预测。
5 结论1) 基于振动理论分析发现，非定常气动力的动态响应过程存在非稳定和稳定两个阶段，传统非定常气动力建模方法的关注点在稳定阶段，即气动力滞环的建模预测，而飞机的真实机动过程常处于进入滞环的非稳定阶段。
2) 采用最小二乘支持向量机进行大迎角非定常气动力建模，通过对某机翼在大迎角范围内做俯仰运动时升力系数、阻力系数和俯仰力矩系数包含初始非稳定过程在内的全过程准确预测，验证了所用方法对非稳定动态过程非定常气动力的预测能力。
3) 通过选择和设计合理的训练输入以激励系统的非线性特征，代替了常见的以多个强迫振动试验数据为输入，简化了训练过程，避免了以滞环预测滞环。在输入的形式选择上考虑到气动力对运动过程的依赖，舍弃了常见的以减缩频率为参数，建立了不依赖于频率的能反映任意运动过程的具有普遍意义的气动力模型。
4) 开始运动的基准状态对气动力的非稳定动态过程存在明显影响。在某些情况下，初始非稳定过程会和稳定滞环重合，即差异表现不明显，但在其他情况下，初始非稳定过程会和稳定滞环存在显著区别。
| [1] | GHOREYSHI M, JIRASEK A, CUMMINGS R M. Computational investigation into the use of response functions for aerodynamic-load modeling[J]. AIAA Journal, 2012, 50(6): 1314-1327. |
| Click to display the text | |
| [2] | MCCRACKEN A J, KENNETT D J, BADCOCK K J, et al. Assessment of tabular models using CFD: AIAA-2013-4978[R]. Reston: AIAA, 2013. |
| [3] | GHOREYSHI M, JIRASEK A, CUMMINGS R M. Reduced order unsteady aerodynamic modeling for stability and control analysis using computational fluid dynamics[J]. Progress in Aerospace Sciences, 2014, 71: 167-217. |
| Click to display the text | |
| [4] | LUCIA D J, BERAN P S, SILVA W A. Reduced-order modeling:New approaches for computational physics[J]. Progress in Aerospace Sciences, 2004, 40(1-2): 51-117. |
| Click to display the text | |
| [5] |
汪清, 钱炜祺, 丁娣. 飞机大迎角非定常气动力建模研究进展[J]. 航空学报, 2016, 37(8): 2331-2347. WANG Q, QIAN W Q, DING D. A review of unsteady aerodynamic modeling of aircrafts at high angles of attack[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(8): 2331-2347. (in Chinese) |
| Cited By in Cnki (6) | Click to display the text | |
| [6] | SILVA W A. Application of nonlinear systems theory to transonic unsteady aerodynamic responses[J]. Journal of Aircraft, 1993, 30(5): 660-668. |
| Click to display the text | |
| [7] |
陈森林, 高正红, 饶丹. 基于多小波的Volterra级数非定常气动力建模方法[J]. 航空学报, 2018, 39(1): 121379. CHEN S L, GAO Z H, RAO D. Unsteady aerodynamics modeling method using Volterra series based on multiwavelet[J]. Acta Aeronautica et Astronautica Sinica, 2018, 39(1): 121379. (in Chinese) |
| Cited By in Cnki | Click to display the text | |
| [8] | GOMAN M, KHRABROV A. State-space representation of aerodynamic characteristics of an aircraft at high angles of attack[J]. Journal of Aircraft, 1994, 31(5): 1109-1115. |
| Click to display the text | |
| [9] | VINOGRADOV Y A, ZHUK A N, KOLINKO K A, et al. Mathematical simulation of dynamic effects of unsteady aerodynamics due to canard flow separation delay[J]. TsAGI Science Journal, 2011, 42(5): 655-668. |
| Click to display the text | |
| [10] |
龚正, 沈宏良. 非定常气动力非线性微分方程建模方法[J]. 航空学报, 2011, 32(1): 83-90. GONG Z, SHEN H L. Unsteady aerodynamic modeling method using nonlinear differential equations[J]. Acta Aeronautica et Astronautica Sinica, 2011, 32(1): 83-90. (in Chinese) |
| Cited By in Cnki (18) | Click to display the text | |
| [11] | WANG Q, HE K F, QIAN W Q, et al. Unsteady aerodynamics modeling for flight dynamics application[J]. Acta Mechanica Sinica, 2012, 28(1): 14-23. |
| Click to display the text | |
| [12] |
汪清, 蔡金狮. 飞机大攻角非定常气动力建模与辨识[J]. 航空学报, 1996, 17(4): 391-398. WANG Q, CAI J S. Unsteady aerodynamic modeling and identification of airplane at high angles of attack[J]. Acta Aeronautica et Astronautica Sinica, 1996, 17(4): 391-398. (in Chinese) |
| Cited By in Cnki (25) | Click to display the text | |
| [13] | GHOREYSHI M, JIRASEK A, CUMMINGS R M. Computational approximation of nonlinear unsteady aerodynamics using an aerodynamic model hierarchy[J]. Aerospace Science and Technology, 2013, 28(1): 133-144. |
| Click to display the text | |
| [14] | BRANDON J M, MORELLI E A. Nonlinear aerodynamic modeling from flight data using advanced piloted maneuvers and fuzzy logic: NASA-TM-2012-217778[R]. Washington, D.C.: NASA, 2012. |
| [15] | WANG Q, QIAN W, HE K. Unsteady aerodynamic modeling at high angles of attack using support vector machines[J]. Chinese Journal of Aeronautics, 2015, 28(3): 659-668. |
| Click to display the text | |
| [16] | CHEN G, ZUO Y, SUN J, et al. Support-vector-machine-based reduced-order model for limit cycle oscillation prediction of nonlinear aeroelastic system[J]. Mathematical Problems in Engineering, 2012. |
| [17] | VAPNIK V. The nature of statistical learning theory[M]. New York: Springer-Verlag, 1995: 161-206. |
| [18] | BURGES C J C. A tutorial on support vector machines for pattern recognition[J]. Data Mining and Knowledge Discovery, 1998, 2(2): 121-167. |
| Click to display the text | |
| [19] | SMOLA A J, SCHOLKOPF B. A tutorial on support vector regression[J]. Statistics and Computing, 2004, 14(3): 199-222. |
| Click to display the text | |
| [20] | SUYKENS J A K, VANDEWALLE J. Least squares support vector machine classifiers[J]. Neural Processing Letters, 1999, 9(3): 293-300. |
| Click to display the text | |
| [21] | MORELLI E, DERRY S, SMITH M. Aerodynamic parameter estimation for the X-43A (hyper-x) from flight data: AIAA-2005-5921[R]. Reston: AIAA, 2005. |
| [22] | MORELLI E. Flight-test experiment design for characterizing stability and control of hypersonic vehicles[J]. Journal of Guidance, Control, and Dynamics, 2009, 32(3): 949-959. |
| Click to display the text | |
| [23] | CHENG C M, PENG Z K, ZHANG W M, et al. Wavelet basis expansion-based Volterra kernel function identification through multilevel excitations[J]. Nonlinear Dynamics, 2014, 76(2): 985-999. |
| Click to display the text | |
| [24] | TRONCHIN L. The emulation of nonlinear time-invariant audio systems with memory by means of Volterra series[J]. Journal of the Audio Engineering Society, 2012, 60(12): 984-996. |
| Click to display the text | |
| [25] |
TISCHLER M B, REMPLE R K.飞机和旋翼机系统辨识: 工程方法和飞行试验案例[M].张怡哲, 左军毅, 译.北京: 航空工业出版社, 2012: 75-76. TISCHLER M B, REMPLE R K. Aircraft and rotorcraft system identification: Engineering methods with flight-test examples[M]. ZHANG Y Z, ZUO J Y, translated. Beijing: Aviation Industry Press, 2012: 75-76(in Chinese). |
| [26] | PRAZENICA R J, KURDILA A J. Multiwavelet constructions and Volterra kernel identification[J]. Nonlinear Dynamics, 2006, 43(3): 277-310. |
| Click to display the text | |
| [27] | CHEN S, GAO Z. Unsteady aerodynamics modeling using Volterra series expansed by basis function[C]//2018 Asia Conference on Mechanical Engineering and Aerospace Engineering. Wuhan: EDP Sciences, 2018. |
| [28] | SCHROEDER M. Synthesis of low-peak-factor signals and binary sequences with low autocorrelation[J]. IEEE Transactions on Information Theory, 1970, 16(1): 85-89. |
| Click to display the text | |
| [29] | IGNATYEV D I, KHRABROV A N. Neural network modeling of unsteady aerodynamic characteristics at high angles of attack[J]. Aerospace Science and Technology, 2015, 41: 106-115. |
| Click to display the text | |
| [30] |
史志伟, 王峥华, 李俊成. 径向基神经网络在非线性非定常气动力建模中的应用研究[J]. 空气动力学学报, 2012, 30(1): 108-112. SHI Z W, WANG Z H, LI J C. The research of RBFNN in modeling of nonlinear unsteady aerodynamics[J]. Acta Aerodynamica Sinica, 2012, 30(1): 108-112. (in Chinese) |
| Cited By in Cnki (22) | Click to display the text | |
| [31] | CAWLEY G C, TALBOT N L C. Preventing over-fitting during model selection via bayesian regularisation of the hyper-parameters[J]. Journal of Machine Learning Research, 2007, 8: 8841-861. |
| Click to display the text | |
| [32] | NASA. CFL3D[CP]. https: //github.com/nasa/cfl3d. |
| [33] | TINOCO E N, BRODERSEN O, KEYE S, et al. Summary of data from the sixth AIAA CFD drag prediction workshop: CRM cases 2 to 5: AIAA-2017-1208[R]. Reston: AIAA, 2017. |
| [34] | MANI M, RIDER B J, SCLAFANI A J, et al. Reynolds-averaged Navier-Stokes technology for transonic drag prediction:A Boeing perspective[J]. Journal of Aircraft, 2014, 51(4): 1118-1134. |
| Click to display the text | |
| [35] | PARK M A, LAFLIN K R, CHAFFIN M S, et al. CFL3D, FUN3D, and NSU3D contributions to the fifth drag prediction workshop[J]. Journal of Aircraft, 2014, 51(4): 1268-1283. |
| Click to display the text |