基于网格速度法的非定常流场模拟和动导数计算

基于网格速度法的非定常流场模拟和动导数计算郭东;徐敏;陈士橹【摘要】Sections 1 through 3 of the full paper explain and evaluate the computation method mentioned in the title, which we believe is effective. Their core consists of; the field velocity or grid velocity approach provides a u-nique methodology for directly calculating aerodynamic responses to step change in flow conditions; the grid time metrics include the velocity caused by the impulsive change in angle of attack but the mesh is not moved accordingly; this approach avoids numerical instabilities and decouples the step change in the angle of attack from a pitch rate; based on this approach, a technique is presented to model longitudinal unsteady flow phenomenon by superposing the step change in the angle of attack upon the impulsive change in pitch rate. In Figs. 4 and 5, numerical results are validated by comparison with experimental results for NACA 0012 airfoil under forced oscillations. To validate further the applicability for the present method, pitch damping derivatives, calculated from the load history of the unsteady flow around a standard research configuration, known as the Basic Finner Missile, are presented in Figs. 10 and 11. Predicted results show indeed good agreement with available wind tunnel data.%基于网格速度法的思想,将迎角的突变和俯仰角速率的突变叠加起来,发展了一套在固定网格上模拟飞行器俯仰振荡非定常流场的方法,该方法不需要实时更新网格,减少了计算时间和所需内存,避免了负体积的出现.首先计算了NACA0006的阵风响应和NACA0012翼型的俯仰振荡,所得结果与实验值和动网格方法符合较好,这表明该方法能够准确模拟此类非定常问题；最后将该方法应用于国外动导数计算标模Basic Finner Missile(BFM)俯仰振荡运动的数值模拟,并计算了其在马赫数1.58 ～2.5的静、动稳定性导数,计算结果与风洞实验值基本吻合,体现了该方法的正确性.【期刊名称】《西北工业大学学报》【年(卷),期】2012(030)005【总页数】5页(P784-788)【关键词】计算流体力学;非定常流场;网格速度法;动导数【作者】郭东;徐敏;陈士橹【作者单位】西北工业大学航天学院,陕西西安710072;西北工业大学航天学院,陕西西安710072;西北工业大学航天学院,陕西西安710072【正文语种】中文【中图分类】V211.3动导数是飞行器动态气动特性设计中的关键参数。

精确导数计算工具是飞行器设计部门所必需的。

随着计算机技术和计算流体力学(CFD)的快速发展，CFD技术已经成为求解飞行器的静、动导数的一种重要手段［1-3］。

非定常流场的准确求解是气动迟滞效应时域仿真的关键，直接影响到静、动导数的计算精度。

目前，利用CFD求解静、动导数时，非定常流场都是采用刚性或者柔性动网格方法来实现的，然而动网格方法需要在每个时间步更新网格，这不仅增加了计算时间和所需内存，还可能出现负体积。

为了求解飞行器迎角或者俯仰角速率的阶跃响应，Baeder等人提出了“网格速度法”，将迎角的阶跃变化，采用对整个流场区域施加一个垂直方向的速度来进行模拟。

在理论上，这只改变了迎角，而对俯仰角速率没有任何影响。

对于角速率的阶跃响应，可以通过对网格点施加一个与旋转中心距离呈线性关系的网格速度来描述，即v(x，y)=˙α·r(x，y)，其中v(x，y)是网格点对应的速度，r(x，y)是网格点相对于旋转中心的距离，从而实现了迎角和俯仰角速率时间历程的解耦，由此得到的结果与精确值吻合较好，而且还不会产生数值振荡［4，5］。

该方法被广泛应用于分析飞行器的阵风响应特性［6-8］。

本文在“网格速度法”的基础上，将迎角的突变和俯仰角速率的突变叠加起来，发展了一套在固定网格上模拟飞行器俯仰振荡非定常流场的方法。

相对“动网格方法”而言，该方法具有以下特点:1)不用实时更新网格，从而减少了计算时间和所需内存;2)不需要考虑网格的密度、类型，避免了负体积的出现;3)可以研究单纯的迎角或者俯仰角速率阶跃响应和阵风响应。

首先通过研究NACA 0006的阵风响应来验证了本文所使用程序的可靠性，然后对NACA 0012翼型俯仰振荡进行了数值研究，通过CT5算例与实验值和动网格方法所得结果进行比较，结果吻合较好，从而表明了该方法和使用的CFD程序能够比较准确模拟此类非定常问题;最后利用该方法数值模拟了国外动导数计算标模Basic Finner Missile(BFM)的非定常俯仰振荡，并计算了其在马赫数1.58˜2.5的静、动导数，计算结果和风洞实验结果吻合较好，体现了方法的正确性。

1 理论基础及数值方法1.1 流场控制方程采用有限体积法离散Euler方程式中，Q=(ρ，ρu，ρv，ρw，ρe)T 为解向量;∂Ω 为流体域Ω的边界;n为边界的外法向矢量，F为对流矢通量。

无粘通量利用Roe通量差分分裂格式离散，采用双时间推进，引入子迭代，利用二阶LUSGS隐式格式推进。

在计算过程中，严格地满足几何守恒定律［9］。

1.2 网格速度法在计算区域内的流场速度可以表示为式中，u、v、w分别是来流速度沿坐标轴的分量，xτ、yτ、zτ是相应的网格速度分量。

当流场中还有外来的速度场(u'，v'，w')时，流场速度变为网格速度方法通过改变网格的速度来模拟这些速度场。

修正后的网格速度可以表示为从而，整个流场的速度为1.3 俯仰振荡模拟方法对于俯仰振荡运动，本文将其视为迎角的突变和角速率的突变两个部分的叠加。

迎角的突变可以用流场在z方向的速度突增来模拟。

此时相当于在流场中注入了一个外部速度场，其3个速度分量分别为式中，α0是初始迎角，Δα是迎角的增量。

对于角速率的变化，可以通过对网格点施加一个与旋转中心距离呈线性关系的网格速度来模拟综上，结合方程(5)和(6)即可得到俯仰振荡运动中网格速度为式中，(xi，yi，zi)是网格点的坐标，(xc，yc，zc)是转轴中心。

基于此，就可以在不使用动网格技术的前提下，模拟飞行器的俯仰振荡运动。

1.4 静、动态稳定性参数辨识对单自由度强迫俯仰振荡，给定简谐振动(无量纲化)形式如下式中，α0为初始迎角，Δα 为振幅，k= ωLref/(2V∞)为减缩频率，t*是无量纲时间，对动态俯仰力矩系数在起始攻角处展开(保留至一阶导数项)如下式中，Cm0为起始攻角处的静态俯仰力矩系数为俯仰静导数为俯仰动导数。

参考文献［3］，利用积分法有2 算例验证为了验证本文的算法和所使用CFD程序的可靠性，首先计算了有实验结果的非定常气动力标准算例。

2.1 NACA 0006翼型迎角的阶跃响应对于迎角的阶跃变化，采用对整个流场区域施加一个垂直方向的速度来进行模拟。

假设翼型初始状态分别为 Ma=0.3，0.5，0.65，0.8，初始迎角为0°，受到垂直阵风的作用，取wg=0.08V∞，导致翼型迎角突然增加约4.6°。

取真实时间步长dt=1.0×10-4s，定义无量纲时间，其中c是参考长度。

图1给出了不同马赫数在此阵风作用下翼型的升力系数随无量纲时间的变化历程。

Lomax［10］推导了线性可压缩流中平板翼型迎角突增α时，阵风响应初始阶段升力系数的理论解图1 迎角阶跃响应下的升力响应图2给出了该理论解与本文计算结果之间的比较，从图中可以看出，在Ma=0.3时，计算结果与理论解吻合得较好，随着马赫数增加，计算结果与理论解的偏差增大，其原因是由于随着马赫数增加，流场的非线性特性增加造成的。

2.2 NACA0012强迫俯仰振荡AGARD CT5是常用的非定常算例，其运动规律为α = α0+ Δαsin(2kt*)，α0=0.016°，Δα =2.51°，k=0.0814翼型绕1/4弦长处做俯仰运动，图3给出了AGARD CT5算例［11］验证所采用的C型网格，网格维数为199×40。

图4给出了升力和力矩系数随时间的变化曲线，图5为不同迎角时翼型上下表面的压强分布，从图中可以发现，计算结果与实验值［12］和动网格方法结果吻合较好，动网格采用刚性旋转法，这表明本文提出的方法和采用的仿真程序在计算此类非定常问题时具有较高的可靠性。

图2 升力响应计算结果与理论解比较图3 CT5算例网格示意图图4 AGARD CT5算例计算结果与实验数据和动网格方法结果比较图5 不同迎角状态物面CP曲线3 BFM静、动导数计算研究Basic Finner Missile(BFM)是国际上用于动导数计算程序验证标准模型［2］，其外形尺寸和计算网格如图6和图7所示。

重心位置为5.0d，Ma=1.58，1.75，1.89，2.10，2.50，减缩频率为 k=0.003 16，初始迎角α0=1.5°，振幅Δα=1.5°。

图8和图9分别给出了Ma=1.58时俯仰力矩系数的响应曲线和迟滞曲线，所得结果与动网格方法吻合较好。