CFD讲义-湍流模型.doc

第三章，湍流模型

第一节， 前言

湍流流动模型很多，但大致可以归纳为以下三类：
第一类是湍流输运系数模型，是Boussinesq于1877年针对二维流动提出的，将速度脉动的二阶关联量表示成平均速度梯度与湍流粘性系数的乘积。即：
3－1
推广到三维问题，若用笛卡儿张量表示，即有：
3－2
模型的任务就是给出计算湍流粘性系数的方法。根据建立模型所需要的微分方程的数目，可以分为零方程模型（代数方程模型），单方程模型和双方程模型。
第二类是抛弃了湍流输运系数的概念，直接建立湍流应力和其它二阶关联量的输运方程。
第三类是大涡模拟。前两类是以湍流的统计结构为基础，对所有涡旋进行统计平均。大涡模拟把湍流分成大尺度湍流和小尺度湍流，通过求解三维经过修正的Navier-Stokes方程，得到大涡旋的运动特性，而对小涡旋运动还采用上述的模型。
实际求解中，选用什么模型要根据具体问题的特点来决定。选择的一般原则是精度要高，应用简单，节省计算时间，同时也具有通用性。
FLUENT提供的湍流模型包括：单方程（Spalart-Allmaras）模型、双方程模型（标准κ-ε模型、重整化群κ-ε模型、可实现(Realizable)κ-ε模型）及雷诺应力模型和大涡模拟。

Zero-Equation Models

One-Equation Models
Spalart-Allmaras
Two-Equation Models
Standard k-e
RNG k-e
Realizable k-e
Reynolds-Stress Model

Large-Eddy Simulation

Direct Numerical Simulation
包含更多
物理机理
每次迭代
计算量增加
FLUENT提供的模型选择
RANS-based
models

湍流模型种类示意图
第二节，平均量输运方程

雷诺平均就是把Navier-Stokes方程中的瞬时变量分解成平均量和脉动量两部分。对于速度，有：
3－3
其中，和分别是平均速度和脉动速度（i=1,2,3）
类似地，对于压力等其它标量，我们也有：
3－4
其中，表示标量，如压力、能量、组分浓度等。
把上面的表达式代入瞬时的连续与动量方程，并取平均（去掉平均速度上的横线），我们可以把连续与动量方程写成如下的笛卡儿坐标系下的张量形式：
3－5
3－6
上面两个方程称为雷诺平均的Navier-Stokes（RANS）方程。他们和瞬时Navier-Stokes方程有相同的形式，只是速度或其它求解变量变成了时间平均量。额外多出来的项是雷诺应力，表示湍流的影响。如果要求解该方程，必须模拟该项以封闭方程。
如果密度是变化的流动过程如燃烧问题，我们可以用法夫雷（Favre）平均。这样才可以求解有密度变化的流动问题。法夫雷平均就是出了压力和密度本身以外，所有变量都用密度加权平均。变量的密度加权平均定义为：
3－7
符号～表示密度加权平均；对应于密度加权平均值的脉动值用表示，即有：。很显然，这种脉动值的简单平均值不为零，但它的密度加权平均值等于零，即：
，

Boussinesq近似与雷诺应力输运模型

为了封闭方程，必须对额外项雷诺应力进行模拟。一个通常的方法是应用Boussinesq假设，认为雷诺应力与平均速度梯度成正比，即：

3－8
Boussinesq假设被用于Spalart-Allmaras单方程模型和双方程模型。Boussinesq近似的好处是与求解湍流粘性系数有关的计算时间比较少，例如在Spalart-Allmaras单方程模型中，只多求解一个表示湍流粘性的输运方程；在双方程模型中，只需多求解湍动能k和耗散率ε两个方程，湍流粘性系数用湍动能k和耗散率ε的函数。Boussinesq假设的缺点是认为湍流粘性系数是各向同性标量，对一些复杂流动该条件并不是严格成立，所以具有其应用限制性。
另外的方法是求解雷诺应力各分量的输运方程。这也需要额外再求解一个标量方程，通常是耗散率ε方程。这就意味着对于二维湍流流动问题，需要多求解4个输运方程，而三维湍流问题需要多求解7个方程，需要比较多的计算时间，对计算机内存也有更高要求。
在许多问题中，Boussinesq近似方法可以得到比较好的结果，并不一定需要花费很多时间来求解雷诺应力各分量的输运方程。但是，如果湍流场各向异性很明显，如强旋流动以及应力驱动的二次流等流动中，求解雷诺应力分量输运方程无疑可以得到更好的结果。

第三节， 湍流模型

3.3.1 单方程（Spalart-Allmaras）模型
Spalart-Allmaras模型的求解变量是，表征出了近壁（粘性影响）区域以外的湍流运动粘性系数。的输运方程为：
3－9
其中，是湍流粘性产生项；是由于壁面阻挡与粘性阻尼引起的湍流粘性的减少；和是常数；ν是分子运动粘性系数。
湍流粘性系数用如下公式计算：

其中，是粘性阻尼函数，定义为：，并且。
湍流粘性产生项，用如下公式模拟：
3－10
其中，，而。其中，和k是常数，d是计算点到壁面的距离；S。定义为：
3－11
由于平均应变率对湍流产生也起到很大作用，FLUENT处理过程中，定义S为：
3－12
其中，，，，平均应变率定义为：
3－13
在涡量超过应变率的计算区域计算出来的涡旋粘性系数变小。这适合涡流靠近涡旋中心的区域，那里只有“单纯”的旋转，湍流受到抑止。包含应变张量的影响更能体现旋转对湍流的影响。忽略了平均应变，估计的涡旋粘性系数产生项偏高。
湍流粘性系数减少项为：
3－14
其中， 3－15
3－16
3－17
其中，，，是常数，。在上式中，包括了平均应变率对S的影响，因而也影响用计算出来的r。
上面的模型常数在FLUENT中默认值为：，，，，，，，。

壁面条件

在壁面，湍流运动粘性设置为零。当计算网格足够细，可以计算层流底层时，壁面切应力用层流应力－应变关系求解，即：
3－18
如果网格粗错不能用来求解层流底层，则假设与壁面近邻的网格质心落在边界层的对数区，则根据壁面法则：
3－19
其中，k=0.419，E=9.793。

对流传热传质模型

在FLUENT中，用雷诺相似湍流输运的概念来模拟热输运过程。给出的能量方程为：
3－20
式中，E是总能量，是偏应力张量，定义为：
3－21
其中，表示粘性加热，耦合求解。如果默认为分开求解，FLUENT不求解处。但是可以通过变化“粘性模型”面板上的湍流普朗特数（Prt），其默认值为0.85。
湍流质量输运与热输运类似，默认的Schmidt数是0.7，该值同样也可以在“粘性模型”面板上调节。
标量的壁面处理与动量壁面处理类似，分别选用合适的壁面法则。
综上所述，Spalart-Allmaras模型是相对简单的单方程模型，只需求解湍流粘性的输运方程，并不需要求解当地剪切层厚度的长度尺度。该模型对于求解有壁面影响流动及有逆压力梯度的边界层问题有很好模拟效果，在透平机械湍流模拟方面也有较好结果。
Spalart-Allmaras模型的初始形式属于对低雷诺数湍流模型，这必须很好解决边界层的粘性影响区求解问题。在FLUENT中，当网格不是很细时，采用壁面函数来解决这一问题。当网格比较粗糙时，网格不满足精确的湍流计算要求，用壁面函数也许是最好的解决方案。另外，该模型中的输运变量在近壁处的梯度要比中的小，这使得该模型对网格粗糙带来数值误差不太敏感。
但是，Spalart-Allmaras模型不能预测均匀各向同性湍流的耗散。并且，单方程模型没有考虑长度尺度的变化，这对一些流动尺度变换比较大的流动问题不太适合。比如，平板射流问题，从有壁面影响流动突然变化到自由剪切流，流场尺度变化明显。

3.3.2 标准模型

标准模型需要求解湍动能及其耗散率方程。湍动能输运方程是通过精确的方程推导得到，但耗散率方程是通过物理推理，数学上模拟相似原形方程得到的。该模型假设流动为完全湍流，分子粘性的影响可以忽略。因此，标准模型只适合完全湍流的流动过程模拟。
标准模型的湍动能k和耗散率ε方程为如下形式：
3－22
3－23
在上述方程中，表示由于平均速度梯度引起的湍动能产生，是用于浮力影响引起的湍动能产生；可压速湍流脉动膨胀对总的耗散率的影响。湍流粘性系数。
在FLUENT中，作为默认值常数，＝1.44，=1.92，，湍动能k与耗散率ε的湍流普朗特数分别为＝1.0，＝1.3。可以通过调节“粘性模型”面板来调节这些常数值。

3.3.3 重整化群κ-ε模型

重整化群κ-ε模型是对瞬时的Navier-Stokes方程用重整化群的数学方法推导出来的模型。模型中的常数与标准κ-ε模型不同，而且方程中也出现了新的函数或者项。其湍动能与耗散率方程与标准κ-ε模型有相似的形式：
3－24
3－25
表示由于平均速度梯度引起的湍动能产生，是用于浮力影响引起的湍动能产生；可压速湍流脉动膨胀对总的耗散率的影响，这些参数与标准κ-ε模型中相同。和分别是湍动能k和耗散率的有效湍流普朗特数的倒数。
湍流粘性系数计算公式为：
3－26
其中，，
对上面方程积分，可以精确得到有效雷诺数（涡旋尺度）对湍流输运的影响，这有助于处理低雷诺数和近壁流动问题的模拟。
对于高雷诺数，上面方程可以给出：，。这个结果非常有意思，和标准κ-ε模型的半经验推导给出的常数非常近似。
在FLUENT中，如果是默认设置，用重整化群κ-ε模型时候是针对的高雷诺数流动问题。如果对低雷诺数问题进行数值模拟，必须进行相应的设置。

重整化群κ-ε模型有旋修正

通常，平均运动有旋时候对湍流有重要影响。FLUENT中重整化群κ-ε模型通过修正湍流粘性系数来考虑了这类影响。
湍流粘性的修正形式为：
3－27
其中，是不考虑有旋计算出来的湍流粘性系数；Ω是FLUENT计算出来的特征旋流数；是旋流常数，不同值表示有旋流动的强度不同。流动可以是强旋或者中等旋度的。FLUENT默认设置＝0.05,针对中等旋度的流动问题，对于强旋流动，可以选择较大的值。
湍动能及其耗散率的有效湍流普朗特数倒数的计算公式为：
3－28
式中，＝1，在高雷诺数流动问题中，，。
湍流耗散率方程右边的R为：
3－29
其中，，，。
为了更清楚体现R对耗散率的影响，我们把耗散率输运方程重写为：
3－30
则： 3－31
在的区域，R的贡献为正；大于。以对数区为例，，，这和标准κ-ε模型中给出的＝1.92接近。因此，对于弱旋和中等旋度的流动问题，重整化群κ-ε模型给出的结果比标准κ-ε模型的结果要大。
重整化群模型中，，。

3.3.4 可实现κ-ε模型

可实现κ-ε模型的湍动能及其耗散率输运方程为：
3－32
3－33
其中，，
在上述方程中，表示由于平均速度梯度引起的湍动能产生，是用于浮力影响引起的湍动能产生；可压速湍流脉动膨胀对总的耗散率的影响。和是常数；，分别是湍动能及其耗散率的湍流普朗特数。在FLUENT中，作为默认值常数，＝1.44，=1.9，＝1.0，＝1.2。
可实现κ-ε模型的湍动能的输运方程与标准κ-ε模型和重整化群κ-ε模型
有相同的形式，只是模型参数不同。但耗散率方程有较大不同。首先耗散率产生项（方程右边第二项）不包含湍动能产生项，现在的形式更能体现能量在谱空间的传输。另外的特色在于耗散率减少项中，不具有奇异性。并不象标准κ-ε模型模型那样把K放在分母上。
该模型适合的流动类型比较广泛，包括有旋均匀剪切流，自由流（射流和混合层），腔道流动和边界层流动。对以上流动过程模拟结果都比标准κ-ε模型的结果好，特别是可实现κ-ε模型对圆口射流和平板射流模拟中，能给出较好的射流扩张角。

湍流粘性系数公式为，这和标准κ-ε模型相同。不同的是，在可实现κ-ε模型中，不再是个常数，而是通过如下公式计算：
3－34
其中，，，，是 is the mean rate-of –rotation tensor viewed in a rotating reference frame with the angular velocity 。模型常数，，而：
，式中W＝，，
我们可以发现，是平均应变率与旋度的函数。在平衡边界层惯性底层，可以得到＝0.09,与标准κ-ε模型中采用底常数一样。

双方程模型中，无论是标准κ-ε模型、重整化群κ-ε模型还是可实现κ-ε模型，三个模型有类似的形式，即都有κ和ε的输运方程，它们的区别在于：1，计算湍流粘性的方法不同；2，控制湍流扩散的湍流Prandtl数不同；3，ε方程中的产生项和Gk关系不同。但都包含了相同的表示由于平均速度梯度引起的湍动能产生，用于浮力影响引起的湍动能产生；可压速湍流脉动膨胀对总的耗散率的影响。
湍动能产生项
3－35
3－36
式中，Prt是能量的湍流普特朗数，对于可实现κ-ε模型，默认设置值为0.85。对于重整化群κ-ε模型，，。热膨胀系数，对于理想气体，浮力引起的湍动能产生项变为：
3－37
在FLUENT程序中，如果有重力作用，并且流场里有密度或者温度的梯度，浮力对湍动能的影响都是存在的。浮力对耗散率的影响不是很清楚，因此，默认设置中，耗散率方程中的浮力影响不被考虑。如果要考虑浮力对耗散率的影响，用“粘性模型”面板来控制。浮力对耗散率影响是用来体现。但并不是常数，而是如下的函数形式：
3－38
v是平行于重力方向的速度分量；u是垂直于重力方向的速度分量。如果流动速度与重力方向相同的剪切流动，＝1，对于流动方向与重力方向垂直的剪切流，＝0。

对于高马赫数的流动问题，可压速性对湍流影响在中体现。

其中，是马赫数，定义为：（是声速）。
默认设置中，只要选择可压速理想气体，可压速效应都是考虑的。
在上述的双方程模型中，对流传热传质模型都是通过雷诺相似湍流动量输运方程得到的。能量方程形式为：
3－39
式中，E是总的能量，是有效导热系数；是偏应力张量，定义为：
3－40
表示的是粘性加热，耦合求解时总是计算。如果不是耦合求解时候，作为默认设置，并不求解该量。如果有需要，需在“粘性模型”面板中设置。
对于重整化群κ-ε模型，有效导热系数为：
3－41
α用（3－28）计算，式中，。事实上，随着的变化而变化，这是重整化群κ-ε模型的一个优点，因为实验中证明，湍流普朗特数随分子普朗特数及湍流而变化。
湍流质量输运处理过程与能量输运过程类似。对于标准κ-ε模型和可实现的κ-ε模型，默认的Schmidt数是0.7，重整化群模型中，是通过方程3－28来计算的，其中，，Sc是分子Schimidt数。

3.3.5 雷诺应力模型（RSM）

雷诺应力模型是求解雷诺应力张量的各个分量的输运方程。具体形式为：

对流项

湍流扩撒项 分子扩散

应力产生项 浮力产生项目

压力应变项 耗散项
3－42
系统旋转产生项

上面方程中，，，，不需要模拟，而，，，需要模拟以封闭方程。下面简单对几个需要模拟项的模拟。
可以用Delay and Harlow [L38]的梯度扩散模型来模拟，即：
3－43
但这个模型会导致数值不稳定，因此FLUENT程序中采用标量湍流扩散模型：
3－44
式中，湍流粘性系数用来计算，根据Lien and Leschziner [L98]，，这和标准κ-ε模型中选取1.0有所不同。
根据Gibson and Launder [L58], Fu [L55], Launder [L88,L89], 压力应变项可以分解为三项，即：
3－45
，和分别是慢速项，快速项和壁面反射项。
，常数。
，，，，。
壁面反射项用于重新分布近壁的雷诺正应力分布，主要是减少垂直于壁面的雷诺正应力，增加平行于壁面的雷诺正应力。该项模拟为：

3－46
式中，＝0.5，，是在垂直于壁面方向上的单位分量，d是到壁面的距离；，，k=0.41。
默认设置时候，FLUENT不计算。如果需要计算时候，在“粘性模型”面板中设置。

线性压力应变模型

对于小雷诺数流动，特别是用双层模型求解近壁流动问题时，FLUENT中通过改进模型常数，，和来改进压力应变项Launder [L91]。这一过程只有在选择双层流模型时候，在“粘性模型”面板上调节。

其中，，参数A和张量不变量A2, A3定义为：

式中，是雷诺应力张量各向异性部分，定义为：
3－47

二阶压力应变模型

二阶压力应变模型由Spezible {L157}等人提出。

3－48

式中，是雷诺应力各向异性张量，定义为：
3－49
平均应变率定义为: ；；模型常数，，，，，，。
二阶压力应变模型不需要考虑壁面反射影响去模拟对数区湍流边界层过程。

浮力对湍流的影响

浮力引起的产生项模拟为：
3－50
其中，Prt是能量的湍流普朗特数，默认设置值为0.85。
对于理想气体，把热膨胀系数的定义代入上式，得：

3－51

35