亚音速下翼型气动性能数值分析

张 猛, 柳兆涛, 姚 程,2, 刘加伟

(1.合肥工业大学 土木与水利工程学院,安徽 合肥 230009; 2.土木工程结构与材料安徽省重点实验室,安徽 合肥 230009)

0 引 言

机翼作为飞机的重要组成部分,机翼翼型的气动性能直接影响飞机的运动特性和使用寿命。不同的翼型有着不同的特点,机翼飞行性能的高低与翼型的选择具有重要的关系,本文通过数值模拟研究翼型的气动性能,可以为一些实际问题提供研究基础,对提升飞机的安全性能和新翼型的研发有一定的意义。在飞行器飞行的过程中,飞行的角度、飞行的速度以及飞行的高度是目前比较关键的因素。由于飞行的高度决定着飞行器周围的大气压,进而影响空气密度,雷诺数会相应地发生改变。对于机翼来说,升力系数和阻力系数是衡量其飞行性能的重要因素。提升飞机的气动性能基本要求就是减小阻力、增加升力和提高升阻比,因此研究机翼最大升力系数、最大升阻比2个变量相当重要。

从20世纪开始国内外研究者便开始对翼型进行一定的研究,研究的内容包括不同的流场因素对翼型气动性能的影响。文献[1]对现有的实验数据进行分析并提取了很多信息,包括对马赫数、雷诺数和攻角的分析,同时更加关注对翼型表面的压力分布和俯仰力矩的分析;文献[2]提出基于低雷诺数翼型的空气动力学的实用研究结果,通过理论分析和实验数据,在合理的范围内解决了理论和实验中出现的各种问题;文献[3]利用Spalart-Allmaras湍流模型对NACA0006翼型和NACA6412翼型进行了气动性能数值研究,分析了这2种翼型在不同的攻角和马赫数的气动性能,研究了攻角和马赫数对翼型升力系数、阻力系数以及升阻比的影响;文献[4]针对3种不同的翼型进行了扰流流动的数值模拟计算,采用SSTk-ω模型得出翼型在固定雷诺数下的升阻系数、升阻比与来流攻角的变化关系,并与翼型的实验数据进行了对比。

本文建立了NACA0012对称翼型的二维模型,研究其在亚音速下翼型气动性能,在亚音速条件下采用湍流模型进行飞行性能模拟并与风洞实验数据对比,分析了雷诺数、马赫数和攻角对翼型的最大升力系数、最大升阻比的影响规律;同时对翼型进行失速特性分析,并在所研究的马赫数范围内给出NACA0012翼型的失速攻角表达式,为翼型失速攻角随马赫数变化的研究提供参考。

1 翼型仿真计算

1.1 湍流模型

不同的湍流模型计算所得到的气动性能是有所差别的,目前在翼型的数值湍流模型的选择上主要为Spalart-Allmaras模型、k-ω模型和SST模型。

根据文献[5]得出SSTk-ω模型更适合用于模拟翼型边界层的流动。该模型为两方程形式,同时集合了k-ε模型和k-ω模型的特点,在近壁面保留了原始k-ω的模型,而在远离壁面的区域采用k-ε模型,其方程式如下:

(1)

(2)

(3)

(4)

(5)

(6)

1.2 计算网格和无关性分析

本文根据美国国家航空咨询委员会(NACA)提供的翼型坐标文件建立翼型的二维模型,通过数值模拟验证发现二维翼型的升阻系数与弦长无直接关系,因此为了方便分析,建立弦长为1 m的数值模型。

同时更好地模拟翼型的流场,本文选择建立C型网格对翼型进行分析,其中网格的边界距离翼型的长度为12.5倍的弦长。

网格的划分决定着计算结果的精度及计算过程的效率,本文在攻角α=0°和马赫数Ma=0.15下进行网格无关性验证[6],验证结果见表1所列。

表1 不同的网格划分比较

通过观察表1中不同网格对升力系数CL和阻力系数CD的影响,最终决定采用505×235的网格。

翼型的计算网格如图1所示。

图1 翼型的计算网格

1.3 边界条件和收敛控制

在实际飞行的状态下,飞机会受到强压缩气流的影响,为提高数值模拟准确性,设置流场的入口边界为压力远场,设置出口边界为压力出口。同时本文采用黏度满足Sutherland定律的理想气体,壁面边界设置为固定壁无滑移条件,环境温度固定为300 K。在二维翼型的模拟计算中,本文采用可压缩流体进行计算,因此选用耦合式求解器。通量采用矢通量分裂Roe格式,动量、湍流动能和湍能耗散率均采用二阶迎风格式。采用湍流强度I和湍流长度尺度l的组合进行模拟,由于湍流强度I的数值与雷诺数密切相关,参考文献[7],计算公式为:

I=0.16Re1/8

(7)

l=0.07L

(8)

其中,L为关联尺寸,对于充分发展的湍流,可取L等于水力直径。

因为本文采用的是基于密度的隐式求解器,所以在进行求解时可选择求解指导器进行迭代运算。Solution Steering可以执行多重网格初始化,不断调整CFL值,保证数值更加快速地收敛。本文在亚音速下且流体可压缩的流动类型下对翼型求解,因此符合Solution Steering对求解类型的要求。

2 模拟结果和分析

2.1 求解模型的准确性验证

NACA0012翼型是目前数值模拟采用较多的低速对称模型,从20世纪开始国外的学者们便对其进行了大量的风洞实验研究[8-9]。

本文选取马赫数为0.15、雷诺数为3.94×106来流条件进行数值模拟,得出不同攻角下的升阻力系数的数值模拟结果,并与文献[10]中的实验数据比较,如图2、图3所示。

图2 升力系数随攻角变化曲线

图3 升力系数-阻力系数曲线

从图2、图3可以看出,翼型模拟得到的升力系数与实验结果有较好的吻合度,而阻力系数与实验数据相比最大误差接近41.5%。通过参考文献[11]可知,阻力系数存在的误差原因有以下几点:

(1) 使用N-S方程计算湍流模型的阻力时结果不精确,尤其是对单元翼型的附着流动。

(2) 单元翼型在大攻角的工况下,对黏性阻力的预测会出现误差。

(3) 对多元翼型来说,采用表面压力和表面摩擦计算的阻力会出现偏差超过100%的情况。

2.2 马赫数和雷诺数对翼型气动性能的影响

本节研究的马赫数为0.10、0.20、0.30、0.40、0.50雷诺数为1×105、5×105、1×106、5×106、10×106,通过数值模拟得出翼型在0°～20°攻角范围内的升力系数和阻力系数,然后将得到的升力系数和阻力系数进行整理用于下文的分析。

2.2.1 马赫数和雷诺数对最大升力系数的影响

升力系数曲线最高点对应的迎角为临界迎角或失速攻角(αs),而对应的升力系数为最大升力系数(CL,max)。从一般意义上来说,衡量一架飞机基本机动性能的优劣,通常看其临界迎角和最大升力系数[12]。最大升力系数是决定飞机起飞着陆性能的重要参数,因此研究翼型最大升力系数在不同马赫数和雷诺数下的变化规律具有重要的意义。

根据绘制的升力系数曲线,由定义取升力系数曲线的最大值为最大升力系数。通过对所有工况的最大升力系数进行归纳,本文得到最大升力系数随马赫数和雷诺数的变化规律如图4所示。

从图4可以看出,在雷诺数为1×105时,翼型最大升力系数(CL,max)为0.8左右,出现在Ma=0.40左右;当雷诺数大于5×105时,同一雷诺数下,翼型最大升力系数变化较大,且最大值出现在低马赫数下;当雷诺数为5×105和1×106时,最大升力系数出现在Ma=0.10左右,为1.17和1.30;当雷诺数为5×106和10×106时,最大升力系数出现在Ma=0.20左右,为1.46和1.59。研究马赫数与最大升力系数的关系得出,在飞机保持飞行高度不变时,提升飞机马赫数会使翼型的失速提前更快,尤其在飞机飞行高度较低时。

2.2.2 马赫数和雷诺数对最大升阻比的影响

升阻比是同一迎角下飞机的升力和阻力的比值。升阻比是衡量飞机空气动力性好坏的重要参数,通过升阻比曲线可以得到最大升阻比Kmax。本文由前面数值模拟得出的升力系数和阻力系数值可得到升阻比并绘制升阻比曲线,曲线可得到升阻比的最大值。研究雷诺数和马赫数对最大升阻比的影响,结果如图5所示。

图5 最大升阻比随马赫数变化曲线

因为升阻比同时考虑了升力和阻力2个因素,所以可以更加全面地反映翼型的飞行性能。从图5可以看出,随着雷诺数的增大,同一马赫数下翼型的最大升阻比增大,这说明在飞机马赫数相同时,随着飞机飞行高度的增加,飞机飞行性能有所下降;在低雷诺数的情况下,翼型的最大升阻比随马赫数增大而先增大后减小,且翼型的最大升阻比出现位置在马赫数为0.20～0.30,这说明在飞机飞行速度不变时,随着飞行高度的降低,飞行器飞行性能在提升。

2.3 翼型最大升力系数和升阻比的攻角位置

根据上节研究的马赫数为0.10、0.20、0.30、0.40、0.50,雷诺数为1×105、5×105、1×106、5×106、10×106情况下,对翼型在0°～ 20°攻角范围内,研究翼型最大升力系数下的攻角α1和最大升阻比下的攻角α2的位置关系,见表2所列,为飞机实际飞行工况提供参考。

表2 翼型最大升力系数和最大升阻比的攻角比较 单位:(°)

从表2可以看出,当雷诺数大于5×105时,在同一马赫数下随着雷诺数增大,翼型最大升力系数和最大升阻比的攻角变化不大;最大升阻比下攻角要比最大升力系数下攻角要低,两者之间的攻角差值为1°～3°,当马赫数为0.10、雷诺数为10×106时,翼型最大升力系数和最大升阻比下的攻角分别为17°和14°;在雷诺数不变的情况下,随着马赫数的增大,翼型最大升力系数和最大升阻比的攻角都在减小。这说明当飞机飞行高度相同时,飞行速度越快,增大翼型攻角,就越容易失速;当飞行速度相同时,飞行高度越大,增大翼型攻角,也容易失速。

2.4 翼型失速特性分析和定量描述

对于机翼的静态失速来说,机翼的迎角达到一定的程度后,翼型上表面的气流会发生流动分离,同时又受外层气流的带动,最后会卷成一个分离涡。因此当机翼的攻角达到一定的程度后就会失速,机翼的升力会减小,此时对应的机翼攻角称为失速攻角。从升力系数曲线上来看,最大升力系数对应的攻角即为失速攻角。

在来流马赫数为0.30和标准大气压下,进行来流攻角为0°～20°的翼型流场数值模拟分析,得到的升力系数如图6所示。

图6 Ma=0.30时升力系数随攻角变化曲线

当攻角为8°时,没有出现翼型分离,这是由于压差阻力很小,绕翼型的阻力绝大部分是由摩擦阻力引起的。当翼型攻角增大为14.01°时,翼型上下翼面压差阻力增大,其尾部出现了分离涡,攻角继续增大,升力系数明显下降,阻力上升,此时翼型进入失速状态。来流攻角为8.00°和14.01°的绕翼型流线分布如图7所示。

图7 来流攻角为8.00°和14.01°的绕翼型的流线分布

本文用曲线拟合来探究翼型在亚音速工况下的失速攻角,为翼型的失速研究提供解决方法。为此较为准确地绘制翼型失速曲线,选取马赫数计算的间隔为0.05,选择翼型失速攻角的工作马赫数为0.15～0.75,静压为标准大气压,求解得出翼型在不同马赫数下的失速攻角与马赫数之间关系,如图8所示。

图8 失速攻角随马赫数变化曲线

曲线拟合的可靠性参考拟合结果报表中的相关系数R2,该数若越接近±1,则表示数据相关度越高,拟合吻合性越好,其表示数据的离散程度,一般要求在0.99以上。通过观察升力系数与失速攻角曲线,可知非线性曲线拟合更加合适,本文通过使用对数函数、指数函数和幂函数进行拟合。从图8可以看出,当马赫数在0.15～0.25范围内,翼型的失速攻角在16°左右;当马赫数大于0.25时,翼型的失速曲线呈现下降趋势,有较强的规律性,可以用曲线进行描述。本文通过3种非线性曲线拟合函数对马赫数为0.25～0.75的翼型失速进行拟合分析,比较R2发现对数函数y=1.64-9.43ln(x-0.02)数据分析框的值更加接近1,拟合得更好,由此可得翼型的失速拟合曲线如图9所示。

图9 翼型失速拟合曲线

通过以上曲线拟合分析,采用分段函数的形式,本节可以定量描述NACA0012翼型在大气压和马赫数在0.15～0.75范围内的失速攻角,即

(9)

3 结 论

本文采用数值模拟研究雷诺数和马赫数对NACA0012翼型的最大升力系数、最大升阻比的影响规律,最后对翼型失速现象进行研究,包括翼型失速特性分析、失速曲线的绘制及定量描述,分析计算结果得到以下结论:

(1) 当Ma=0.15、Re=3.94×106时,模拟分析获得的翼型升力系数与实验数据十分吻合,而阻力系数有一定的误差,因此分析了阻力系数存在偏差的原因。

(2) 当雷诺数为5×106、10×106时,最大升力系数随马赫数的变化波动较大,且变化趋势基本相同,最大升力系数出现在Ma=0.20左右分别为1.46和1.59,为所研究范围飞机的最佳飞行状态;在低雷诺数的情况下,翼型的最大升阻比随马赫数增大而先增大后减小。且翼型的最大升阻比出现位置在马赫数为0.20～0.30,为飞机在高空飞行提供飞行参数依据;翼型最大升阻比下攻角比最大升力系数下攻角低,两者之间的攻角差值为 1°～3°。

(3) 通过对翼型在亚音速下的失速特性进行分析,提出了翼型的失速攻角随马赫数变化的函数表达式,通过定量描述法,可以确定翼型在某一马赫数范围内的失速攻角值,为研究亚音速下飞行器的安全飞行范围提供了参考。