航空航天
关键郭正
( 国防科学技术大学航天科学与工程学院,长沙410073)
摘要层流转捩到湍流及翼型表面分离泡的产生是绕翼型低雷诺数流动的重要特征,对流场产生很大影响。针对这一流动现象,选取E387翼型为研究对象,采用求解雷诺平均N-S 方程的有限体积法,结合当前应用较广泛的S p a l a r t-A ll m a r as、S S T k-ω、Rea li za b l e k-ε三种湍流模型,对雷诺数为3 ×105 时翼型绕流流场进行了数值模拟;并将结果与Se li g等人的风洞实验数据进行对比,评估三种湍流模型对绕翼型低雷诺数流动的模拟能力。基于对翼型阻力计算不准原因的分析,提出了一种基于M i c h e l转捩判据的数值模拟方法,分别从分离泡的模拟、壁面摩阻分布、翼型阻力系数等方面与实验数据进行对比。结果表明该方法可以较好地模拟低雷诺翼型绕流流场。 关键词低雷诺数分离泡转捩湍流模型翼型
V211.3;
中图法分类号文献标志码 A
翼型在低雷诺和不可压来流时的气动性能越来越得到人们的关注。这是由于近年来许多领域的研究热点都与低雷诺数流动有着紧密的关系。在人们对长期驻留高空的长航时无人机开展研究时,如何提高其机翼在低速时的气动特性,以及如何改善发动机叶片在高空时的工作效率都成为人们关心的问题,究其本质都是研究翼型在低雷诺数流动下的气动特性的问题。对低雷诺数下翼型绕流流场进行仿真分析时,层流转捩过渡到湍流是绕翼型低雷诺数流动的一个重要特征;这一流动现象对气动特性有很大影响,在一些工况下会产生层流分离、湍流再附、分离泡等一系列复杂的流动现象,这些现象都对数值模拟提出了挑战。
当前在计算流体力学( C o m p ut a t i o n a l Fluid Dy- n a m i cs,CFD)领域中Re yn o l ds平均方法( RA N S)以其可以有效地避免大量的计算,更经济可行地模拟湍流而得到广泛应用;但该方法并不能很好地处理转捩问题,如何对这一问题进行建模国内外许多学者都对此做了研究。R. B L a n g t r y 等人对两方程湍流模型SS T k-ω进行了改进,提出了一种可以预测转捩的两方程湍流模型[1];并对S809 翼型进行了仿真计算。B ac h a等人也基于Spa l a r t-A ll m a r as湍流模型,结合史密斯-加姆波尼的e9关联曲线的转捩点判据,提出了一种可以用于转捩预测的湍流模型,并针对NLF-0416 翼型进行了算例计算[2]。国内学者扬中也于2010 年提出了一种基于k-ω模型添加控制层流动能的三方程湍流模型对低压涡轮叶栅
[3]
进行了算例验证。
针对上述问题,本文拟在评估现有湍流模型的
基础上通过考虑转捩以提高对低雷诺数流动的模
拟效果。选取雷诺数在10 量级E387( E)翼型绕流
5
流场进行数值研究。这一研究对象及条件的选取主要基于以下两个原因: E387 (E)自20 世纪60 年
代被设计出来以后在滑翔机上取得了成功,后来逐
渐成为一种典型的低雷诺数风洞测试翼型,先后在
德国的Stuttg art 及NASA 的L a n g l e y 和U I UC 均做过翼型吹风试验,有大量试验数据可以对比;雷诺数
5
在10 量级时的流动包含前文所述的较复杂的流动
现象因而具有一般性。
文中选用的湍流模型为工程中得到广泛应用的Spa l a r t-A ll m a r as,SS T k-ω、Rea li zab l e k-ε湍流模型,结合RA N S方法数值模拟流场,得到攻角在-6°
2013 年5 月3 日收到
第一作者简介:关键( 1988—),男,硕士。研究方向:计算流体力学。E-ma i l:go d o l ph i n@s i na.c o m。
智能家居系统
7276
科 学 技 术 与 工 程 13 卷
至 12°时升阻 力 系 数随攻角变化曲线,以 此 作 为 依 据评估三种湍流模型对低雷诺 流场仿真的模拟能 力,并对结 果 进 行 分 析,进而探索一种基于 M i c h e l 转捩判 据的湍流数值模 拟 方 法,并从分离泡的模 拟、壁面摩阻分布、翼型阻力系 数 等 方 面 与 Se li g 的 风洞实验数据进行对比,验证其可行性。
理定律,在标准 k-ε 湍流模型基础提出了改进模型 即 Rea li zab l e k-ε 湍流模型,其中 k 和 ε 的输运方程 如下
[5]
( ρk ) +
( ρku i ) = [ ( + μt ) k ] + G
- ρε;
μ k
t
x i x j x
j
σ k
( ρε) +
( ρεu i ) =
[ (
+
) ]
+ ρC E ε - μt ε
μ 1 t
x i x j
x j
σ ε
1 控制方程与湍流模型
ε2
ρC 2
。
k + 槡v ε
1. 1 主控方程
流体在运动时必须遵守流动控制方程,根据物 质的连续性、动量守恒得到流体的控制方程
光标跟随
1. 4 SST k-ω 湍流模型
SST k-ω 湍 流 模 型 是 在 Menter V l co x 提 出 的 标 准 k-ω 湍流模型的基础上,结合 k-ε 湍流 模 型 的 优 点建立的一个分区模型; 充分利用 k-ω 模型对逆压 梯度比较敏 感,能 够模拟较大分离运动的特点,在 靠近壁面的附面层采用 k-ω 湍流模型; 在远离附面 层的流场中,采用 k-ε 湍流模型,克服了 k-ω 模型对 自由来流比较敏感的缺点,提高了模型稳定性。该 模型能够模拟中等分离流动,尤其在逆压梯度较大 的流场模 拟 中 表 现 出 。该两方程湍流模型的控 制方程为
[6]
:
- ρ + "·( ρv) 连续方程: = 0。
t
动量方程:
-
-
- - - -
= - "·(
) + ρg + F 。
( ) + "·( ) t ρv ρv ·v
τ -
式中,τ 为应力张量,
定义如下式 -
τ = μ [
( "·v + "·v
T )
- -
]
2 -
-
·v I 。 " 3 -
式中,v 为分子黏性,I 为单位张量。 u i D ρk ρωk + [ ( μ + σk μt ) x ]
;
k 1. 2 Spa l a r t -A ll ma r a s 湍流模型
S -A 是一 种 求 解动力涡黏输运方程的简单模 * = τij x + β Dt D ρω x j j j u i - βρω2 + γ [ ( μ + σ μ ) ω ]
+ 型,它包含了 一 组 新发展的单方程模型,方 程 里 不 去计算局 部 剪 切 层厚度相关的长度尺度。该 模 型 是专门用于求解航空领域的壁面限制流动,对于受 逆压梯 度作用的边界层 流 动,已取得比较好的效 果,原始的 S-A 模型是一种低雷诺数模型,由于引入 壁面函数的方法在粗壁面条件 下也可以取得很好 的应用。该湍流模型的控制方程为
[4]
= i j Dt v t x j x j ω t x j
1 k ω
2ρ( 1 - F 1 ) σω
。
ω x x j j
2 低雷诺下全湍流模型数值仿真及分析
2. 1 计算对象及网格
为方便与 Se li g 等人的实验数据进行比对,数值 -
D v 5
仿真选取 实验时的雷诺数即 3 × 10 ,计 算 对 象 为 E387( E ) 翼 型,使 用 P o i nt w i se 前 处 理 软 件 对 E387 ( E ) 翼型进行“C ”型拓扑结构网格划分。为提高计 算效率 又不损失计算精度,对网格进行无关性分 析,针对上 述 C 型结构网格分别对网格单元数为 3 × 104 ~ 9 × 104 的 7 套 网 格 进 行 计 算,获 取 攻 角 1. 5°时的升力系数并对结果分析,由 图 1 所 示 升 力 系数随网格单元数的变化曲线中可以看出,当网格 单元数目为 7 × 104 时升力系数随网格量的增加不 在有显著变化,因此本文计算域最终确定网格单元
- - = c b 1 ( 1 - f t 2 ) Sv + 1
- - "v + [ σ "·( v + v)
Dt 2 (
c w1 f w - 2 t )
c - c b2 ( "v) b 1 2
2 - + f t1 ΔU 。 ]
2
κ
1. 3 Re a li zab l e k-ε 湍流模型
与一方程模型不同,在湍动能 k 的方程基础上, 再引入一关于湍动耗散率 ε 的方程便形成了两方程 模型,称为标准 k-ε。该模型是由 Launder 和 Spa l d i n g 提出的。但由于标准 k-ε 对时均应变率特别大的流 动可能导致负的正应力,为使流动更符合湍流的物
24 期关键,等:绕翼型低雷诺数流动的数值仿真7277
凯膜过滤技术数为7 ×104,高度为20 倍弦长,第一层网格高度保
证y + = 1 以满足不同湍流模型对网格的要求,计算
域及壁面附近网格见图 2 和图3。
边界条件采用速度入口边界及压力出口边界条
件。利用S I MP L EC 算法求解流场,采用二阶精度空
间离散格式,湍流模型分别选用S-A、SST k-ω、R
数字振镜ea li zab l e k-ε进行计算。
图3 壁面附近网格
下,S-A 模型本身就具有模拟低雷诺数流动的能力,
而本文所选两方程模型均为高Re 数模型,计算结果
差于S-A 模型。
当攻角大于8°时,两方程模型的升力估计均优
于S-A 模型,这是由于S-A 这种单方程模型对于长
度尺度变化不敏感,尤其在分离时壁面约束流动转
变为类自由剪切流时计算结果有偏差。
阻力预测方面,通过图 5 分析三种湍流模型计
算结果,在低雷诺数来流绕翼型流动问题的仿真
中,S-A 湍流模型阻力预测结果要优于其他两种湍
流模型,但总的来说三种湍流模型对于阻力的估计
均不如对升力估计的准确,这是完全湍流模型一种
固有缺陷,即完全湍流模型无法考虑层流附面层向
湍流附面层的转捩过程,无法捕捉到层流分离转捩
至湍流及湍流再附和翼型表面分离泡等流动现象,
致使在层流区也采用湍流模型计算,势必阻力估计
不准,因此如果想提高仿真的准确性,必须考虑转
捩过程。
图1 升力系数随网格单元数变化曲线
图2 E387( E)翼型C 型网格
2. 2 结果与分析
采用三种湍流模型数值计算了E387 ( E)翼型
在- 6°~12°的升阻力系数,用以评估这三种湍流
模型对低雷诺问题的模拟能力。
通过图 4 可以看出在升力的计算方面,攻角小
于8°时三种湍流模型对于升力的预测与实验数据
吻合较好,其中以S-A 模型对当前雷诺数下升力估
计较准,这主要源于当壁面附件网格较密的情况
考虑转捩的仿真方法建立
3
如果要提高数值模拟精度,就必须考虑低雷诺
层流向湍流的转捩过程,通过加入转捩判据,采用
将层湍流分开计算方法,建立一种可以判转捩的数
值模拟方法,下面就分别从其转捩的判断和层流湍
流的处理两方面对该方法进行阐述。
3. 1 转捩点预测
7278
科 学 技 术 与 工 程 13 卷
合 求得曲面层流动量厚 度,基 于 大 量 地 观 测 Th wa i t es 发现动量厚度 θ 与无黏势流速度 u e ( x) 有 如下关系:
u e ( x ) d u e ( x ) 。 2
d ( θ2 ) = A - B θ υ dx υ dx
式中,A 、B 根据 Th wa i t es 推荐取 A = 0. 45,B = 6 为最 佳经验值。
对上式积分得到边界层动量厚度
x
0. 45υ ∫
θ2
u 5
( x) dx 。 =
6
e u e ( x) x = 0
基于此法得到的 θ( x) 对于绝大多数的不可压 层流流动条件下,其误差度可以控制在 ± 3% 内。
3. 2 层流湍流的处理
由于转捩点后湍流发展区比较短,工程上普遍 认为转捩点后即可认为流动已 经为完全发展的湍 流,在预测的 转 捩 点处将流动区域分开两部分,在 转捩点以前区域采用层流模型计算,转捩点以后区 域采用湍流模型计算,以此来模拟层流转捩为湍流 的过程,其中湍流模型选用前文仿真结果较好的 S -A 模型,利用此 方 法 分别从分离泡现象的仿真,壁 面 摩阻的分布及阻力系数的影响 三方面进行算例计 算,评估方法可行性。
图 4 三种湍流模型计算得到升力系数曲线
与实验数据对比
4 算例计算及结果分析
4. 1 分离泡现象模拟
图 5 三种湍流模型计算得到阻力系数曲线
与实验数据对比
本文转捩判断基于 M i c h e l 方法和 G a m bo r n i 的 e 9 关联曲线的结合。由 S i mth 给出的转捩时的 Re θt r [7]
与 Re x 之间有如下关系 :
= 1. 174 (
1 + )
Re tr
22 400 0. 46 6
Re θt r ; 0. 1 × 10 ≤x Re x tr
6 Re x ≤ 40 × 10 。
式中,Re x 为当地雷诺数,Re θ 为基于动量厚度的雷
图 6 Re = 3 × 105
,5°攻角时 E387( E ) 上
表面油流显示图像
诺数,上式适用于 0. 1 × 106
≤ Re ≤
40 × 106 x 从 Se li g 的油流显示实验可以看出,在 Re = 3 ×
攻 角 为 5° 时,E387 ( E ) 翼 型 上 表 面 的 流 动 在
不可压绕翼 型 流 动,其中来流速度给定后,可 以 通 5
过求解 Eu l e r 方 程 得 到 无黏流动速度分布 u e ( x ) 10 。 其中动量厚度可由 Th wa i t es 积 分 法 求 得,该 方 法 适
40% 弦 长左右处层流发生分离,在 图 6 中 40% ~
24 期
关 键,等: 绕翼型低雷诺数流动的数值仿真
7279
油管锚
55% 弦长处的斑点状图像正是实验前油雾在翼型上 表面形成的图像,说明此区域流动近乎滞止即属于 分离区已发 生 分 离。在 55% 弦长处有一条很明显 的分界线,这 是 油 滴在回流作用下形成的,这 一 条 线位于再附点之前,
58% 弦长处出现一条不明显的 分界线,这是 由 于 再 附 后,表面油滴被向上游的部 分回流及向下游的流体冲刷所至,因而此处是流动 再附点,在再附线下游就湍流流动区域。
针对上述从实验中观察到的流动现象,本文分 别采用层流、Spa l a r t-A ll m a r as 湍流模型、及所建立的 考虑转捩的湍流数值模拟方法对流场进行仿真,即
在 Re = 3 × 105
,
5° 攻 角 时 对 E387 ( E ) 翼 型 进 行 计 算,并将计算得到的压力云图及翼型壁面附近流线 图给出。
图 7 所示即为全湍流流动仿真结果,壁面附近 流线平滑,无 涡 结 构且流动并未发生分离,与 油 流 显示实验给出的结果不符,没有捕捉到真实流动中 分离泡现象,而图 8 所示全层流计算结果显示,流动 在后缘处已 经 发 生 分 离,形成了不稳定的涡结构, 流动过程中始终伴随着涡的生成和脱落,这与真实 流场有很大的区别。
图 8 Re = 3 × 105
,
5°攻角 E387( E ) 层流流场 图 9 Re = 3 × 105
,
5°攻角 E387( E ) 固定点
转捩法仿真得到流场
图 7 Re = 3 × 105
,5°攻角 E387( E ) 湍流流场
烘干炉图 10 翼型表面分离泡
考虑了转捩 的 方 法 计 算得到流场见图 9,以 整 体上来看流场是稳定的,这与图 7 计算所得流场近 似相同,但在 机 翼 上表面处有一个小分离泡,分 离 泡比较薄具体细节如图 10 所示,不同于全层流模型 计算时在尾缘处发生分离有很多不稳定的涡,用这
种方法计算得到的结果显示分 离后的流动很快就 发生了再附,可 以 捕捉到流动的分离泡结构的,比 单纯使用湍流模型计算和层流 计算得到的结果更 符合真实流动。