三种典型翼型边界层稳定性对比分析

(
f' f － f″ ξ ξ
)
（ 2）
f = f' = 0 ； η = η e ， f' = 1 。 其中导数表 边界条件为 η = 0 ， m 表示无量纲梯度： m = 示对 η 求导， εm x + + ， b = 1 + εm 。 ε = ， L m v 1． 2 线性稳定性方程 扰动的发展模式可以分成时间模式和空间模式， 时 间模式的扰动在空间某点仅随时间规律变化， 空间模式 的扰动传播则沿空间某一方向变化， 而后者是可以通过 实验结果验证的
， 为之后的相
［4 ］
关研究提供了丰富的数据 。 朱万琳等人在二维不可压 流动的稳定性数值分析基础上进行了转捩的判定 。 为了便于工程上翼型的选取及优化， 针对多个典型翼型 边界层稳定性的对比分析具有重要的参考意义 。 首先 本文采用萘升华实验方法验证对称翼型 SD8020 数值模 拟的准确性， 并拟定稳定性计算的幅值放大因子的值， 然后针对常规翼型 （ NACA0012 ） ， 自 然 层 流 翼 型 （ NACA64204 ） 和超临界翼型 （ RAE2822 ） 这三种典型翼型 采用线性理论方法计算在不同弦长雷诺数情况下转捩 位置的变化， 研究流向压力梯度对转捩点雷诺数 、 转捩 并分析影响 点弦长百分数和其相应的扰动频率的影响， 这三种翼型稳定性增长的原因 。
S 方程） 计算扰动波在流 利用线性稳定性方程（ O场中的传播， 由于在二维基本流条件下三维扰动的失稳 点大于二维扰动的失稳点， 因此引入的二维扰动也是合 理的
［7 ］
捩预测 。
， S方 即假定 β = 0 。 这里二维流无量纲化的 O［1 ］
程可以表示为一个四阶常微分方程 － 2 α ″ + α =
［3 ］
变换基本流方程
［5 ］
， 取η=（
u e / vx ） y， 这里 u e 表示自 槡 （ 1）
y） 的无量纲流函数为： 由流速度 。相对于流函数 ψ（ x， f（ η） = y） ψ（ x， u vx 槡e m +1 f f″ + m［ 1 － （ f' ） 2］ 2 = ξ f'
由此可得连续方程和动量方程为： （ bf″） ' +
482
空
气
动
力
学
学
报
N
第 29 卷
y ） = （ y ） e i （ αx + β z － ωt ） ψ（ x，
（ 3）
对应上述实验条件， 采用 e 方法对 SD8020 翼型进 行转捩预测 。 图 2 为翼型上表面扰动增长曲线， 计算 20 条不同频率的扰动增长曲线， 选取常用的扰动放大 因子 N = 9 。当某一频率增长率最先达到 N = 9 时认为 用虚线与频率增长率曲线的交点对应 计算得到转捩点， 的弦长位置来表示 。 计算得到转捩位置在 79% 左右，
N
0
引
言
［1 ］
1
1． 1 ， 因此
［2 ］
数值模拟方法
边界层平均流求解 S法 采用无粘流和边界层方法求解基本流， 利用 F-
摩擦阻力可占亚音速运输机总阻力的 50% ， 而层 流的阻力比相同雷诺数湍流的阻力要小 90% 减小表面摩擦阻力可以有效地降低飞行器的总阻力 。
为了计算机翼的摩擦阻力从而进一步更合理地预测机 翼飞行特性， 转捩位置的预测成为了重要前提工作 。 对 于高空飞行外界扰动较小的特点， 线性稳定性方法的 e N 方法仍然是目 前 工 程 应 用 中 最 常 用 的 转 捩 预 测 方 Falvy 通过修正 Mack 的 法。L． J． Runyan 和 D． George稳定性方法计算了翼型在不同工况下的转捩位置并与 得到相应的最大增长率的值 实验对照，
第 29 卷 第 4 期 2011 年 8 月
空 气 动 力 学 学 报 ACTA AERODYNAMICA SINICA
Vol． 29 ， No． 4 Aug． ， 2011

1825 （ 2011 ） 04048105 文章编号： 0258-
图5 Fig． 5
6 三种翼型在弦长雷诺数 Re c = 9 × 10
下不同频率 N 因子曲线及转捩位置 N factor curves and transition locations of three airfoils for Re c = 9 × 10 6
484
空
气
动
力
学
学
报
第 29 卷
不稳定区域也变大 。相比高弦长雷诺数， 低弦长雷诺数 下相同频率的扰动更早进入中性曲线的不稳定区域， 即 且沿频率增长方向扰动不稳定区范围 扰动增长的更早， 更小 。因此随弦长雷诺数的增长， 同一翼型最不稳定扰 动波频率也变大 。
三种典型翼型边界层稳定性对比分析
王 菲， 额日其太， 王 强， 苏沛然
（ 北京航空航天大学 能源与动力工程学院，北京 100191 ）
摘
N 要： 为了便于工程上翼型的选取及优化， 本文针对三个典型翼型边界层稳定性进行对比分析 。 首先采用 e 方
法对翼型 SD8020 进行转捩预测， 并采用萘升华实验法检验数值计算的准确性 。再利用该数值方法对比了三种典型 不同弦长雷诺数下的压力系数分布 、 扰动增长率以及最不稳定扰动波的频率 。 分析结果表 翼型上表面的转捩位置、 NACA0012 最先发生转捩， 204 和 RAE2822 都保持着较长的层流区； 压力梯度 明， 在相同弦长雷诺数下， 而 NACA64转捩点雷诺数变大， 而不是工程上常采用的固定值 。 对扰动增长有很大影响； 同一翼型随弦长雷诺数增长， 关键词： 线性稳定性； e 方法； 转捩； 翼型； 压力梯度 中图分类号： V211． 4 文献标识码： A
9 与实验结果吻合， 因此以下方法均采用 e 方法进行转
这里 表示流函数的型函数 。 α 和 β 分别表示流向和 ω 表示扰动的圆频率 。 对应于空间模 展向的扰动波数， 式α， β 为复数， α = α r + iα i ， β = β r + iβ i ； ω 为实数 。 这 里定义无量纲扰动相速度 c： c= ω = c r + ic i α （ 4）
6 9 诺数（ Re c = 9 × 10 ） 下扰动增长曲线， 采用 e 方法进行
转捩预测 。三个翼型中 NACA0012 的转捩点雷诺数最 也即在相同弦长雷诺数下最早发生转捩 。由于在二 小， S 波对转捩的影响起主导作用， 维流动中 T而顺压梯度
图1 Fig． 1 萘升华实验
S 波有稳定作用 。 对应 于 图 4 压 力 系 数 分 布， 又对 TNACA0012 仅有 10% 弦长的顺压区， 远小于后两者， 因
6 表 1 统计了三种翼型在弦长雷诺数分别为 6 × 10 、
9 × 10 6 、 12 × 10 6 下的转捩点雷诺数和对应弦长百分数 。 从表中可以看出， 工程上常采用雷诺平均方法计算全湍 流流场
［8 ］
或固定转捩点雷诺数来假设转捩位置的方法
Baidu Nhomakorabea
是不准确的 。同一翼型随着弦长雷诺数的增长， 转捩点 雷诺数也随之变大， 但并未按线性比例增长， 其对应的 弦长百分数反而减小， 也就是说转捩点相对提前 。 由于 转捩位置的变化会对翼型性能计算产生很大的影响 区是十分必要的 。
3
计算结果及分析
N 利用上述 e 方法， 选用三种典型翼型 NACA0012 ，
NACA64204 ，RAE2822 （ 图 3 ） ， 计算来流马赫数 Ma = 0． 3 ， 攻角为 0 的条件下这三种翼型压力分布及转捩位 置。图 4 为三种翼型上表面压力系数 C p 分布曲线， 其 204 和 RAE2822 分别对应自然层流翼型和 中 NACA64超临界翼型， 这两种翼型压力分布特点是都具有较长的 在前缘附近剧烈的顺压梯度之后， 有较长 顺压梯度区， 204 和 一段弱顺压区， 最后达到最低压力点 。 NACA64RAE2822 最低压力点分 别 在 弦 长 的 40% 和 55% 位 置 处。 图 5 （ a） （ b ） （ c ） 为三种翼型上表面在相同弦长雷
［7 ］
。
204 翼 图 6 （ a） 、图 5 （ b ） 和图 6 （ b ） 表示 NACA646 6 6 型分别在弦长雷诺数为 6 × 10 、9 × 10 和 12 × 10 不
同频率扰动增长率 。 由这三幅图可以看出随弦长雷诺 数增长， 在层流区内由于顺压区与逆压区的长度比变 其转捩点雷诺数也变大 。 当弦长雷诺数增大到 12 大， × 10 6 时， 在顺压区扰动增长率已经达到了 N = 9 ， 因此 可以判定在最小压力点之前转捩已经发生 。 NACA64204 翼型扰动增长的特点是扰动通过翼型前缘后部大 范围弱顺压梯度后， 逆压梯度对其影响非常明显， 即使 很小的逆压梯度也会使低频扰动迅速增长 。 弦长雷诺数的增长使高频扰动增长率变大， 因此使 最不稳定频率（ 对应转捩位置的扰动频率） 也相应提高 （ 图 7） ， 而不同翼型相同弦长雷诺数下对应最不稳定频 率也存在很大差别 。
［6 ］
x du e u ， f' = ， ξ= u e dx ue
， 所以本文采用更接近实验条件的空
间模式方法进行扰动增长计算 。将流函数定义为：
*
收稿日期： 2010 － 03 － 11 ；
修订日期： 2010 － 06 － 01
mail： wangfei@ sjp． buaa． edu． cn 作者简介： 王菲（ 1981 － ） ， 男， 辽宁沈阳人， 博士生， 主要从事流动稳定性数值计算方面研究 ． E-
2
实验验证
本实验在北京航空航天大学流体力学研究所低速
将萘溶于丙酮 风洞中完成 。实验利用萘的易升华特性， 溶液并喷涂于机翼模型表面 。 由于转捩点前后气流对 壁面的剪切力变化很大， 使萘在湍流区内升华速度加 因此模型表面湍流区的白色萘涂层先消失， 从而在 快， 层流和湍流之间形成了明显的边界 。 选择实验状态为 Re c = 7． 8 × 10 5 ， Ma = 0． 087 ， 攻角为 － 2° 以保持较长的 该翼型转捩位置在 78% 弦长左右（ 图 1 ） 。 层流区，
iv 2 4
： （ 5）
iαR（ u － c） （ ″ － α2 ） － iαRu″ ue c ， 这里 c 表示弦长 。 v
边界条件为 y = 0 ， = 0， ' = 0 ； y = δ， = 0， ' = 0 。R 表 Re c = 示弦长雷诺数 Re c ，
N 转捩预测采用工程上最为常用的 e 方法预测翼型
［9 ］
，
可见准确计算翼型的转捩点并分别求解层流区和湍流
图8 Fig． 8
NACA64204 不同弦长雷诺数下的中性曲线 The neutral curves for different Reynolds numbers for NACA64204
4
结
论
本文 首 先 采 用 萘 升 华 实 验 方 法 验 证 对 称 翼 型 SD8020 数值模拟的准确性， 并拟定 eN 方法的幅值放大 然后采用该数值方法分别计算了在不同的弦 因子的值， 长雷诺数下三种典型翼型的压力分布， 扰动增长率以及 计算了翼型上表面转捩位置， 并 最不稳定扰动波频率， 分析了不同翼型压力分布对扰动发展和转捩位置的影 响。结果表明： （ 1 ） 在相应实验条件下采用幅值放大因 子 N = 9 的转捩预测位置与实验结果吻合 。 （ 2 ） 在三种 NA相 同 弦 长 雷 诺 数 下， 典型翼 型 转 捩 预 测 结 果 中， CA0012 最先发生转捩， 204 和 RAE2822 都 而 NACA64保持着较长的层流区 。 （ 3 ） 压力梯度对扰动增长有很 如图 8 所示， 同一翼型， 随弦长雷诺数增长， 中性曲 线位置向雷诺数更大频率更高的方向偏移， 中性曲线的 大影响， 顺压梯度对二维翼型扰动有明显的稳定作用 。 （ 4 ） 同一翼型随弦长雷诺数增长， 顺压区内扰动增长更
上表面的转捩位置 。扰动增长率 N 为： N=
∫
x
x0
－ α i dx
（ 6）
图2
5 翼型 SD8020 在弦长雷诺数 Re c = 7． 8 × 10 下
其中 － α i 为给定频率下的空间放大率 。
不同频率 N 因子曲线及转捩位置 Fig． 2 N factor curves and transition locations of SD8020 airfoils for Re c = 7． 8 × 10 5
Naphthalene sublimation experiment
第4 期
王
菲等： 三种典型翼型边界层稳定性对比分析
483
此其扰动也越不稳定 。 三种翼型的高频扰动增长起始 位置总是比低频更靠近前缘， 而受顺压梯度影响， 高频 增长速度也更容易减慢甚至衰减 。 相对的低频扰动增 长点起始位置更靠后， 而其受顺压梯度影响的区域也减 少。若低频扰动增长率在进入逆压区时并未开始衰减， 则逆压将会导致扰动迅速增长 。 这种现象与中性曲线 的分布规律是一致的 。 正是由于逆压梯度使速度剖面 产生拐点， 有拐点的速度剖面的稳定性界限比无拐点的 剖面低得多（ 拐点准则）