笛卡尔网格在气动设计中的应用研究

罗剑波，栗莉，方明恩，罗帅，綦龙，张岩

（航空工业洪都，江西南昌330024）

笛卡尔网格在气动设计中的应用研究

罗剑波，栗莉，方明恩，罗帅，綦龙，张岩

（航空工业洪都，江西南昌330024）

采用非贴体笛卡尔网格及有限体积法求解Euler方程的方法对两个模型的气动特性进行了计算，并考虑粘性项对阻力结果进行了修正，将修正后的计算结果与风洞试验值进行对比分析，结果表明：修正后的CFD计算值与风洞试验值吻合较好，验证了在气动设计中笛卡尔网格方法具有计算精度高、求解速度快的特点，适用于气动方案初步设计阶段的实际工程应用。

笛卡尔网格；Euler方程；数值模拟；粘性阻力修正

0 引言

在飞行器气动设计过程中，尤其在气动方案初步设计阶段，往往需要进行大量CFD仿真，以获得各状态下飞行器的气动特性，用于气动方案的选型。以往采用的方法是：首先对设计的气动外形生成结构/非结构网格，再求解N-S方程或Euler方程，最后对求解结果进行处理分析。采用这种方法通常能获得较为精确的解，但网格划分和流场求解过程需花费大量的人力和计算时间，限制了其在气动方案设计阶段的应用。

近年来，基于笛卡尔网格技术求解Euler方程的方法越发成熟，与传统CFD仿真方法不同的是：该方法首先在全流场域生成各向尺寸一致的粗糙网格，再根据模型结构在物面附近自动逐步加密，得到尺寸合适的流场网格，最后求解Euler方程得到流场结果，该方法网格生成效率高、流场求解速度快，能大大缩短方案阶段的气动外形设计时间。

1 Cart3D程序简介

本文采用的笛卡尔网格划分及求解Euler方程的过程在Cart3D程序中完成。Cart3D程序能通过定义网格区域及网格密度，自动捕捉模型的几何特征，快速生成非贴体笛卡尔网格（图1），极大地压缩网格生成时间；求解模块基于有限体积法求解Euler方程，采用显式龙格库塔法时间推进，空间离散采用迎风格式，并采用通量限制器提高精度，具有TVD高阶格式特性，在求解过程中，可采用多重网格、自适应网格等措施提高计算效率和计算精度，适用于亚、跨、超音速复杂飞行器的外流场分析。

2 控制方程与数值算法

三维可压缩流动的Euler方程可以表示为：

式中：U为守恒变量；F、G、H为无粘通量，分别表示为：

式中：ρ、p、e分别表示密度、压力和单位质量总能；u、v、w分别为x、y、z三个方向的速度分量。假设来流为完全气体，则方程组可以使用下式封闭：

对于空气，γ取1.4。

3 算例验证

3.1 ONERA M6机翼

ONERA M6机翼因其外形简单，在跨音速环境中机翼表面绕流呈现出局部超音速流动、激波、边界层分离等复杂的流动状态，因而成为验证CFD方法的典型算例。对ONERA M6机翼划分的笛卡尔网格，如图2所示，采用Euler方程求解，并将计算的压力分布与NASA标准计算结果及ONERA S2MA的风洞试验值进行对比分析，验证采用笛卡尔网格离散并求解Euler方程方法的计算精度。计算条件为：马赫数Ma=0.84，攻角α=3.06°，侧滑角β=0°，Re=11.72×106（基于平均气动弦长bA=0.64607m），参考面积Sref=0.7532m2。

图3左侧是采用笛卡尔网格求解Euler方程得到的ONERA M6机翼上表面压力云图，右侧是NASA求解N-S方程得到的标准计算结果，可以看出两种计算方法得到的上翼面压力分布具有相当高的一致性，压力分布结果基本吻合，采用笛卡尔网格的方法很好地捕捉到了弹翼上表面的λ型激波和波后的流动。

提取ONERA M6机翼展向z/B=0.20、0.44、0.65、0.80、0.90、0.95处6个典型截面位置的压力分布，并与ONERA S2MA风洞的试验值进行对比。如图4所示，图例中Exp表示风洞试验值，Euler为采用笛卡尔网格方法的计算值。结果表明：弹翼下表面流动相对简单，笛卡尔网格方法精确地模拟了弹翼下表面的流动，各截面处计算的压力值与试验的压力值吻合较好；而弹翼上表面由于出现激波、边界层分离等复杂流动现象，笛卡尔网格方法的计算结果与试验值之间有一定的偏差，弹翼上表面前缘处的压力峰值的位置及压力峰值的模拟较为准确，但对弹翼上表面中部激波位置及激波强度的模拟不大准确；总体上看计算结果与试验值吻合较好，可验证笛卡尔网格方法的计算精度。

3.2 某飞行器

采用笛卡尔网格方法对某飞行器算例的升阻特性随攻角和马赫数的变化趋势进行大量数值计算，并将计算与风洞试验结果进行比较。飞行器头部附近网格如图5所示。

利用单机计算了某飞行器在亚音速、超音速下共24个飞行状态的升阻特性，完成整个网格划分及流场求解过程总计耗时不超过3H。由于笛卡尔网格方法仅求解Euler方程，不能预测物面的摩擦阻力，故在计算完成后对阻力系数计算结果进行了粘性修正。

本文对阻力系数采用的粘性修正方法为：

1）利用Fluent等求解N-S方程，得到飞行器的零升阻力系数Cd0NS；

2）利用笛卡尔网格方法求解Euler方程，得到零升阻力系数Cd0E；

3）计算零升阻力系数差量ΔCd0=Cd0NS-Cd0E；

4）将零升阻力系数差量ΔCd0叠加到笛卡尔网格方法的阻力系数求解结果，得到修正后的阻力系数。

图6～图8为飞行器升阻力系数的计算值与试验值的对比，其中：图例Exp表示风洞试验值，Euler为笛卡尔网格方法的CFD计算值，CFD阻力修正为对笛卡尔网格方法的CFD阻力计算值进行粘性修正的结果。可以看出：升力系数的CFD计算值与试验值基本重合；而未修正的阻力系数计算值与试验值存在差量，但粘性修正后的阻力系数曲线与风洞试验值基本重合，且对飞行器的阻力发散马赫数的预测较为准确，说明阻力修正方法可行。由此可见，笛卡尔网格方法能较为准确的预测飞行器的升力系数，阻力系数也能通过粘性修正获得较为精度的结果。

与传统采用贴体结构/非结构网格进行空间离散，并采用Fluent、CFX等求解器求解N-S方程的方法比较，采用笛卡尔网格求解Euler方程的方法不仅具有较高的计算精度，而且因其网格生成迅速、求解效率高，能大幅度缩短设计周期。就上述飞行器模型而言，利用单机对24个飞行状态从建模、求解到后处理，采用笛卡尔网格求解Euler方程的方法共耗时不超过3H，而相同条件下若利用传统方法处理耗时则不少于4天。

4 结论

通过对ONERA M6机翼及某飞行器模型的计算验证可知：

1)笛卡尔网格离散方法能较为准确地模拟物面的流动，包括激波位置、激波强度及物面流动分离现象，具有较高的精度，可满足方案阶段气动特性快速设计的工程需求；

2)利用笛卡尔网格离散方法求解Euler方程，配合多重算法的加速方法，网格生成迅速、求解速度快，计算结果可靠，即使在单机情况下，也比常规采用N-S方程求解器的效率提高至少10倍，从而大大缩短设计周期。

[1]肖涵山，陈作斌，刘刚，江雄。基于Euler方程的三维自适应笛卡尔网格应用研究[J]。空气动力学学报，2003；21（2）：202-210.

[2]戚姝妮，张祖庚，董军。基于自适应笛卡尔网络之三维标模气动特性的数值模拟[J]。气动研究与试验，2005；22（3）：1-8.

Application Research of Cartesian Grid in Aerodynamic Design

Luo Jianbo,Li Li,Fang Mingen,Luo Shuai,Qi Long,Zhang Yan
(AVIC-HONGDU,Nanchang,Jiangxi 330024)

This paper adopts the way of non-body-fitted Cartesian grid and finite volume method for solving Euler equation to calculate the aerodynamic characteristics,and corrects drag result with giving consideration to viscosity term.The corrected calculation result is analyzed by comparing with wind tunnel test value.The result shows:the corrected CFD calculation value fits wind tunnel test value well,so it demonstrates that Cartesian grid method features highcalculation accuracy and fast solving speed,and is applicable to the practical engineering application of aerodynamic scheme in initial design stage.

Cartesian grid；Euler equation；numerical simulation；correction of viscous drag

2017-08-29）

>>>作者简介 罗剑波，1989年5月出生，2011年毕业于南京航空航天大学，工程师，现主要从事飞行器气动设计工作。