飞行器动导数是指动稳定性导数^[1]，它是飞行器绕重心转动所引起的气动力系数和力矩系数对相对运动参数的导数，是飞行器导航系统，控制系统以及动态品质分析不可缺少的原始气动参数，直接影响着飞行器的飞行品质，特别是飞行器机动飞行或操纵时的飞行品质。目前,动导数获取手段有理论分析计算、风洞实验和飞行试验。飞行试验和风洞实验费用大、周期长。理论分析计算包括工程估算以及数值模拟。以一些经验以及半经验公式为代表的工程估算方法是早期研究动导数的主要途径，而这些方法的局限性以及计算机技术的快速发展，加上在非定常流动数值模拟方面的进展使得结合N-S/Euler方程的数值计算成为研究动导数计算的热点和趋势，涉及的理论包括牛顿－玻尔兹曼理论，内伏牛顿理论，基于爆炸波比拟和激波－膨胀波比拟的杂交方法，活塞理论，非定常升力面理论等^{[2, 3, 4, 5, 6, 7]}。国外，Erdal Oktay发展了非定常欧拉方法^[8]，Soo Hyung Park建立了基于双时间推进的DADAI 方法^[9]，Scott M Murman用非线性折算频率方法^[10]模拟单一频率的强迫振动引起的响应求解动导数；国内的研究工作也取得了显著的成果，袁先旭^[11]、孙涛^[12]研究了减缩频率对动导数的影响，史爱明^[13]、卢学成^[14]也是用非结构动态网格重构技术计算了超声速导弹的动导数，而孙智伟^[15]、将胜炬^[16]使用了定常的N-S方程对动导数的计算方法开展了研究。综合来看，这些工作对动导数的数值求解起到很大的推动作用，但是由于计算的复杂性，方法不是十分系统化，精度也较低，因而动导数的数值计算还需要做进一步大量研究。

本文结合CFD计算技术发展了飞行器超声速动导数计算方法。以国外动导数标模^[17]Finner导弹为例，基于滑移网格技术，对动态气动设计特性中的滚转阻尼导数，以及俯仰组合动导数进行了模拟，并采用多参考系方法，对滚转阻尼导数进行了数值计算。计算结果表明，无论是基于滑移网格的非定常方法，还是结合多参考系模型的准定常方法，都具有比较高的计算精度。 1 计算方法 1.1 控制方程

文中采用的控制方程为三维积分形式的雷诺平均N-S方程。

点击浏览公式

空间离散采用二阶迎风格式——通量差分分裂(Roe-FDS)格式，采用双时间步推进并用隐式LU-SGS格式迭代求解，应用当地时间步长、残值光顺、预处理和多重网格加速收敛。湍流模型采用对逆压梯度流动模拟精度较高的k-ω SST湍流模型。 1.2 滑移网格

滑移网格是一种特殊的动网格，它将问题的求解域划分成多个计算域，即动区域和静止区域。各计算域之间采用定义好的分界面进行关联，动区域沿着滑移面进行运动。在滑移网格问题中，动区域运动是相对于静止参考系进行跟踪的，因此，没有运动参考系附加在计算域上，简化了穿过分界面的通量传递。与一般动网格相比，网格的运动仅沿着滑移面进行整体滑移，没有网格变形重构的过程，相对节省了产生新网格所需的计算量，而且运动时各区域网格质量不发生变化，可以提高计算精度。 1.3 多参考系模型

多参考系模型(MRF)是不同旋转或移动速度的每个区域的稳态近似。同样它将计算域分成多个小的子域，每个子域可以有自己的运动方式，或静止或旋转或平移。此方法适用于边界上流动区域几乎均匀混合的情况。以同时存在静止和旋转域情况为例，在两个子域上分别建立静止坐标系和旋转坐标系，控制方程在每个子域分别求解，在交接面上通过将速度换算成绝对速度形式进行流场信息交换。由于进行的是稳态近似，多参考系模型是一种很高效的方法。 1.4 动导数计算方法 1.4.1 强迫振动法

强迫飞行器做小幅度俯仰振荡时，依据傅立叶公式展开，刚体飞行器所受非定常气动力矩可以写为表达式：

点击浏览公式

依据气动力导数的概念，上述小幅度强迫俯仰振动的刚体飞行器所受非定常气动力矩还可以表示为：

点击浏览公式

当刚体飞行器以低频ω做小幅度振荡时，其模型运动方程可以简化为：

点击浏览公式

点击浏览公式

点击浏览公式

假定飞行器以两种不同角速度ω_z1,ω_z2等速上仰运动至相同攻角，依据飞行力学小扰动理论对其力矩进行展开并略去高阶量可得：

点击浏览公式

点击浏览公式

点击浏览公式

前面提出的强迫振动法以及差分法都是非定常方法，可以用来求解纵向，横航向动导数。而一般的飞行器，具有顺流方向的面对称性或者线对称性，横向气动力在0°迎角，定速转动状态下一般保持不变，且横向更多的是考虑稳定滚转特性，因而通常使用准定常方法进行滚转动导数的计算。

飞行器进行小角速度等速滚转时，令驱动力矩为M_x，阻尼力矩为M_kx，根据小扰动理论：

点击浏览公式

点击浏览公式

点击浏览公式

本文使用的基于滑移网格以及多参考系模型的求解方法要求将流场划分为动域以及静域，两个域之间通过交接面进行数据传递。动域中的网格刚化后沿着交接面进行滑移运动，以此来实现物面的运动。本文涉及到的边界条件有远场边界，物面边界以及交接面边界。其中远场边界使用基于Riemann不变量的无反射边界条件，物面边界无滑移，绝热，交接面使用基于通量守恒的交接面边界。 2 动导数计算公式无量纲化

一般来讲，计算动导数的公式都是无量纲的，因此，需要对本文的几种方法进行无量纲化处理。 2.1 强迫振动法

计算模型为国际标准动导数模型Finner导弹，图 1为计算用的Finner导弹模型及网格，其中d=1m，采用结构网格，根据本文的方法将流场分为动域和静域，两部分通过interface交接面连接，这个交接面的网格不要求完全一致，给复杂模型网格划分带来方便。最终划分的网格总量为390万，其中，动域网格量260万，静域130万。

图 1 计算模型及网格Fig. 1 Computational model and grid

图选项

使用非定常方法求解动导数的关键在于准确地预测瞬时力矩。本文选取M_∞=1.58，α₀=0°为计算状态，减缩频率为k=0.0158226，振幅为α_m=1°，也即振荡规律可写为

图 2为计算得到的瞬时力矩随迎角的变化曲线。可以明显地看出，本文的非定常强迫振动方法能够很好地预测迎角与力矩系数之间的迟滞效应，计算的俯仰力矩系数迟滞环与文献中的非常一致。

图 2 力矩系数随迎角变化规律图Fig. 2 Time history of pitching moment coefficient with angle of attack in pitching cycle

图选项

根据前面得到的动导数无量纲计算公式，可求解俯仰组合动导数，计算结果见表 1。

文献中的俯仰组合动导数计算值为-506，误差为3.62%。由此可见本文的计算方法对俯仰组合动导数的计算是比较准确的。

一般使用强迫振动方法计算动导数多使用积分法进行结果辨识，本文则使用稳定后的平衡位置处单点气动数据计算动导数。表 2是两种方法的对比，可以看出，两者均能比较准确的得到动导数结果，相比下，本文的单点数据计算量小，能够方便快速得到结果。

为了进一步就本文的小幅度强迫振动方法进行验证，分别计算在马赫数为1.58，1.74，1.88，2.1，2.55五个状态下的俯仰组合动导数。

图 3的结果表明，重心位置在X_g=5.0d时，本文所计算的动导数与风洞试验数据较为一致，反映出在超音速范围内，随着马赫数的增加，俯仰组合动导数的逐渐减小，且趋势变缓和。由于目前的技术限制，动导数试验结果的误差带相当大，与试验结果相比，本文采用数值模拟技术较为准确地计算了Finner导弹的俯仰组合动导数，在文中给定的马赫数范围内很好地模拟出动导数的变化趋势。

图 3 不同马赫数下的俯仰组合动导数Fig. 3 Pitching dynamic stability derivative versus Mach number

图选项

本文的计算初始迎角均为0°，导弹以5°/s,10°/s的角速度分别等速上仰至5°。然后根据得到的瞬时力矩即可求出俯仰组合动导数。计算结果见表 3。

表 3的结果表明，使用差分法也能比较准确地计算俯仰组合动导数。

给定导弹以5°/s的速度绕轴滚转，使用前面所述的准定常方法，计算得到的滚转阻尼导数为:

由表 4可以看出，准定常虽然计算量小，耗时短，但是精度较低。而强迫振动以及差分方法精度比较高。对于动导数的数值计算，关键在于计算非定常气动力，非定常方法基于这一点提出，因而精度也稍高，相对而言，准定常方法是一种近似的方法，对流场特性细节的捕捉能力稍弱，所以精度较低。

采用基于滑移网格的非定常方法以及基于多参考系模型的准定常方法对Finner导弹标准模型动导数进行了计算研究和分析。相比以往的使用网格重构进行计算的方法，滑移网格的引入使得计算量减小，并且由于不存在网格变形带来的影响，因此可以提高求解的精度。由本文的算例可以看出所发展的非定常方法以及准定常方法能够较为准确地计算俯仰及滚转动导数。其中小幅度强迫振动法以及差分法两种非定常方法的关键在于准确预测瞬时力矩；准定常方法虽然精度稍差，但是计算量小，快捷。同时，本文的三种方法不仅可以用于超音速动导数辨识，也同样适用于亚音速以及跨音速的动导数计算。

本文采用的动导数辨识方法给出了与风洞试验较为吻合的数值计算结果，以及一致的动导数的变化趋势，表明所采用的动导数计算方法的正确性，具备较好的工程应用前景。