LES,DNS,RANS三种模拟模型计算量比较及其原因

LES,DNS,RANS模型计算量比较
摘要：湍流流动是一种非常复杂的流动，数值模拟是研究湍流的主要手段，现有的湍流数值模拟的方法有三种：直接数值模拟（Direct Numerical Simulation: DNS），Reynolds平均方法（Reynolds Average Navier-Stokes: RANS）和大涡模拟(Large Eddy Simulation: LES)。

直接数值模拟目前只限于较小Re数的湍流，其结果可以用来探索湍流的一些基本物理机理。

RANS方程通过对Navier-Stokes方程进行系综平均得到描述湍流平均量的方程；LES方法通过对Navier-Stokes方程进行低通滤波得到描述湍流大尺度运动的方程，RANS和LES方法的计算量远小于DNS，目前的计算能力均可实现。

关键词：湍流；直接数值模拟；大涡模拟；雷诺平均模型
1 引言
湍流是空间上不规则和时间上无秩序的一种非线性的流体运动，这种运动表现出非常复杂的流动状态，是流体力学中有名的难题，其
性。

传统计算复杂性主要表现在湍流流动的随机性、有旋性、统计[]1
流体力学中描述湍流的基础是Navier-Stokes（N-S）方程，根据N-S 方程中对湍流处理尺度的不同，湍流数值模拟方法主要分为三种：直接数值模拟（DNS）、雷诺平均方法（RANS）和大涡模拟（LES）。

直接数值模拟可以获得湍流场的精确信息，是研究湍流机理的有效手段，但现有的计算资源往往难以满足对高雷诺数流动模拟的需要，从而限制了它的应用范围。

雷诺平均方法可以计算高雷诺数的复杂流动，但给出的是平均运动结果，不能反映流场紊动的细节信息。

大涡模拟基于湍动能传输机制，直接计算大尺度涡的运动，小尺度涡运动对大尺度涡的影响则通过建立模型体现出来，既可以得到较雷诺平均方法更多的诸如大尺度涡结构和性质等的动态信息，又比直接数值模拟节省计算量，从而得到了越来越广泛的发展和应用。

2 直接数值模拟(DNS)
湍流直接数值模拟(DNS)就是不用任何湍流模型，直接求解完整的三维非定常的N - S 方程组,计算包括脉动在内的湍流所有瞬时运动量在三维流场中的时间演变。

2.1控制方程
用非稳态的N - S 方程对紊流进行直接计算, 控制方程以张量形式给出: 21i i i j i j i j i u u u p u f v t x x x x ρ∂∂∂∂+=-+∂∂∂∂∂ （1）
=0i j
u x ∂∂ （2） 2.2主要数值方法
由于最小尺度的涡在时间与空间上都变化很快，为能模拟湍流中的小尺度结构，具有非常高精度的数值方法是必不可少的。

2.2.1谱方法或伪谱方法
所谓谱方法或伪谱方法是目前直接数值模拟用得最多的方法，简单来说，就是将所有未知函数在空间上用特征函数展开，成为以下形式：
()()()()123
(,)mnp m n p m n p
V x t a t x x x ψϕχ=∑∑∑ （3） 其中m ψ，n ϕ与p χ，都是已知的正交完备的特征函数族。

在具有周期性或统计均匀性的空间方向一般都采用Fourier 级数展开，这是精度与效率最高的特征函数族。

在其它情形，较多选用Chebyshev 多项式展开，它实质上是在非均匀网格上的Fourier 展开。

此外，也有用Legendre, Jacobi, Hermite 或Laguerre 等函数展开，但它们无快速变换算法可用。

如将上述展开式代入N-S 方程组，就得到一组()mnp a t 所满足的常微分方程组，对时间的微分可用通常的有限差分法求解。

在用谱方法计算非线性项例如V ϖ⨯的Fourier 系数时，常用伪谱法代替直接求卷积。

伪谱法实质上是谱方法与配置法的结合,具体做法是先将两量用Fourier 反变换回到物理空间，再在物理空间离散的配置点上计算两量的乘积，最后又通过离散Fourier 变换回到谱空间。

在有了快速Fourier 变换(FFT)算法以后，伪谱法的计算速度高于直接求两Fourier 级数
的卷积。

但出现的新间题是存在所谓“混淆误差”，即在做两个量的卷积计算时会将本应落在截断范围以外的高波数分量混进来，引起数值误差。

严重时可使整个计算不正确甚至不稳定，但在多数情形下并不严重，且有一些标准的办法可用来减少混淆误差，但这将使计算工作量增[]2加。

2.2.2高阶有限差分法
高阶有限差分法的基本思想是利用离散点上函数值
i f 的线性组合来逼近离散点上的导数值。

设i F 为函数()j f x ∂∂的差分逼近式,则
j j j
F f α=∑ （4） 式中系数j α 由差分逼近式的精度确定，将导数的逼近式代入控制流动的N - S 方程,就得到流动数值模拟的差分方程。

差分离散方程必须满足相容性和稳定性。

2.3 优点
(1)直接数值求解N-S 方程组，不需要任何湍流模型，因此不包含任何人为假设或经验常数。

(2)由于直接对N - S 方程模拟,故不存在封闭性问题,原则上可以求解所有湍流问题。

(3)能提供每一瞬时三维流场内任何物理量（如速度和压力）的时间和空间演变过程，其中包括许多迄今还无法用实验测量的量。

(4)采用数量巨大的计算网格和高精度流体力学计算方法，完全模拟湍流流场中从最大尺度到最小尺度的流动结构，描写湍流中各种尺度的涡结构的时间演变，辅以计算机图形显示，可获得湍流结构的清晰与生动的流动显示。

2.4 缺点
DNS 的主要缺点是要求用非常大的计算机内存容量与机时耗费。

据Kim ,Moin &Moser 研[]3
究 ,即使模拟Re 仅为3300 的槽流,所用的网点数N 就约达到了6210⨯ ,在向量计算机上进行了250 h 。

3 雷诺平均模拟(RANS)
雷诺平均模拟（RANS ）即应用湍流统计理论, 将非稳态的N - S 方程对时间作平均,求解工程中需要的时均[]4
量。

所谓湍流模式理论，就是依据湍流的理论知识、实验数据或直接数值模拟结果，对Reynolds 应力做出各种假设，即假设各种经验的和半经验的本构关系，从而使湍流的平均Reynolds 方程封闭。

3. 1 控制方程
对非稳态的N - S 方程作时间演算, 并采用Boussinesp 假设,得到Reynolds 方程
''
21i j i i i j i j j j j j u u u u u p u f v t x x x x x ρ∂∂∂∂∂+=-+-∂∂∂∂∂∂ （5）
=0i i u x ∂∂ （6）
式中,附加应力可记为''ij i j pu u τ=-,并称为雷诺应力。

这种方法只计算大尺度平均流动,而所有湍流脉动对平均流动的影响,体现到雷诺应力ij τ中。

正因为雷诺应力在控制方程中的出现,造成了方程不封闭，为使方程组封闭,必须建立模型。

3. 2 主要方法
目前工程计算中常用的湍流模型从对模式处理的出发点不同，可以将湍流模式理论分类成两大类：一类引入二阶脉动项的控制方程而形成二阶矩封闭模型，或称为雷诺应力模型，另一类是基于Boussinesq 的涡粘性假设的涡粘性封闭模式，如零方程模型，一方程模型和二方程模型。

3.2.1雷诺应力模型
雷诺应力模型(RSM)从Reynolds 应力满足的方程出发,直接建立以''i j u u 为因变量的偏
微分方程, 将方程右端未知的项（生成项，扩散项，耗散项等）用平均流动的物理量和湍流的特征尺度表示出来，并通过模化封闭。

封闭目标是雷诺应力输运方程:
''
''''''''i j
i j j
i k i k j k ij ij ij k k k u u u u u u u u u u u D t x x x φε∂∂∂∂+=--++-∂∂∂∂ （7） 式中ij φ 是雷诺应力再分配项, ij D 是雷诺应力扩散项, ij ε是雷诺应力耗散[]5项。

典型的平均流动的变量是平均速度和平均温度的空间导数。

这种模式理论，由于保留了Reynolds 应力所满足的方程，如果模拟的好，可以较好地反映Reynolds 应力随空间和时间的变化规律，因而可以较好地反映湍流运动规律。

因此，二阶矩模式是一种较高级的模式，但是，由于保留了Reynolds 应力的方程，加上平均运动的方程整个方程组总计15个方程，是一个庞大的方程组，应用这样一个庞大的方程组来解决实际工程问题，计算量很大，这就极大地限制了二阶矩模式在工程问题中的应用。

2.2.2涡粘性模型
在工程湍流问题中得到广泛应用的模式是涡粘性模式。

这是由Boussinesq 仿照分子粘性的思路提出的，即设Reynolds 应力为， ,,,22()33i j T i j j i k k ij ij u u U U U k νδδ=-+++
（8） 这里12
i j k u u =是湍动能，T ν称为涡粘性系数，这是最早提出的基准涡粘性模式，即假设雷诺应力与平均速度应变率成线性关系，当平均速度应变率确定后，六个雷诺应力只需要通过确定一个涡粘性系数T ν就可完全确定，且涡粘性系数各向同性，可以通过附加的湍流量来模化，比如湍动能k ，耗散率ε，比耗散率w 以及其它湍流量ετ/k =，ε/2/3k l =，k q =，根据引入的湍流量的不同，可以得到不同的涡粘性模式，比如常见的ε-k ，k-w 模式，以及后来不断得到发展的τ-k ，q-w ，k-l 等模式，涡粘性系数可以分别表示为
ενμ/2k C T = ，ωνμk C T =，τνμk C T =，ωνμ2q C T =，l k C T μν=（9）
3. 3 优点
(1) 对计算机的要求较低，同时可以得到符合工程要求的计算结果。

(2)一旦给定合理的Reynolds 应力模型，可以很容易地从RANS 方程解出湍流的统计量，所需要的计算资源小。

(3)几乎能对所有雷诺数范围的工程问题求解，并得出一些有用的结果。

3. 3 缺点
(1)对不同类型的湍流，需要采用不同的Reynolds应力模型，甚至对于同一类型的问题，对应于不同的边界条件需要修改模型的常数。

(2) 由于不区分旋涡的大小和方向性，对旋涡的运动学和动力学问题考虑不足，不能用来对流体流动的机理进行描述。

(3) 对于非定常流动、大分离流动、逆压力梯度数值模拟等问题，受湍流模型条件的限制，很难得到满意的计算结果。

(4)严重依赖流场形状和边界条件，普适性差，计算很大程度上依赖于经验。

4 大涡数值模拟(LES)
湍流大涡数值模拟（LES）是有别于直接数值模拟和雷诺平均模式的一种数值模拟手段。

利用次网格尺度模型模拟小尺度紊流运动对大尺度紊流运动的影响即直接数值模拟大尺度紊流运动, 将N-S方程在一个小空间域内进行平均（或称之为滤波），以使从流场中去掉小尺度涡，导出大涡所满足的方程。

4. 1 基本思想
湍流运动是由许多大小不同的旋涡组成的。

那些大旋涡对于平均流动有比较明显的影响，而那些小旋涡通过非线性作用对大尺度运动产生影响。

大量的质量、热量、动量、能量交换是通过大涡实现的，而小涡的作用表现为耗散。

流场的形状，阻碍物的存在，对大旋涡有比较大的影响，使它具有更明显的各向异性。

小旋涡则不然，它们有更多的共性，更接近各向同性，因而较易于建立有普遍意义的模型。

基于上述物理基础，LES把包括脉动运动在内的湍流瞬时运动量通过某种滤波方法分解成大尺度运动和小尺度运动两部分。

大尺度要通过数值求解运动微分方程直接计算出来,小尺度运动对大尺度运动的影响将在运动方程中表现为类似于雷诺应力一样的应力项,该应力称为亚格子雷诺应力,它们将通过建立模型来模拟。

实现大涡数值模拟,首先要把小尺度脉动过滤掉,然后再导出大尺度运动的控制方程和小尺度运动的封闭方程。

4. 2 过滤函数
大涡模拟第一步就是把一切流动变量划分成大尺度量和小尺度量，这一过程
称之为滤波。

滤波运算相当于在一定区间内按一定条件对函数进行加权平均，其目的是滤掉高波数而只保留低波数，截断波数的最大波长由滤波函数的特征尺度决定。

目前较为常用的滤波函数主要有以下三种：Deardorff 的盒式(BOX)滤波函数、富氏截断滤波函数和高斯(Gauss)滤波函数。

不可压常粘性系数的紊流运动控制方程为N-S 方[]6
程: j
ij i j j i i x S x P x u u t u ∂⋅∂+∂∂-=∂∂+∂∂)2(1γρ （10） 式中：S 拉伸率张量，表达式为：2/)//(i j j i ij x u x u S ∂∂+∂∂=；γ分子粘性系数；ρ流体密度。

设将变量i u 分解为方程(11)中i u 和次网格变量(模化变量)'i u ,即'+=i i i u u u ,i u 可以采用leonard 提出的算式表示为:
（11）
式中)(x x G '-称为过滤函数，显然G(x)满足
4. 3 控制方程
将过滤函数作用与N-S 方程的各项，得到过滤后的紊流控制方程组：
(12) 由于无法同时求解出变量i u 和j i u u ，所以将j i u u 分解成i j i j ij u u u u τ=⋅+，ij τ即称为次网格剪切应力张量(亦称为亚格子应力)。

由此动量方程又可写成：
(13)
式中ij τ代表了小涡对大涡的影响。

x d x u x x G x u i i '
''-=⎰+∞
∞-)()()(⎰+∞
∞-=1)(dx x G ()(2)1i j ij i j i j u u S u P t x x x γρ∂∂⋅∂∂+=-+∂∂∂∂()(2)1i j ij ij i j i j j u u S u P t x x x x τγρ∂⋅∂∂∂∂+=-+-∂∂∂∂∂
4. 4 常用亚格子模式及其特点
目前，在大涡模拟中经常广泛采用的亚格子模型有标准的Smagorinsky 模型、动态涡粘性模型、动态混合模型、尺度相似模型、梯度模型、选择函数模型[]7等。

其中Smagorinsky 模型被广泛应用。

4.4.1亚格子涡粘和涡扩散模型
不可压缩湍流的亚格子涡粘和涡扩散模型采用分子粘性和分子热扩散形式，即
kk ij ij t ij S τδντ3
12+= (14) i
t i x T ∂∂=θκ (15) 以上公式中t ν和t κ分别称为亚格子涡粘系数和亚格子涡扩散系数；(1/2)[(/)(/)]ij i j j i S u x u x =⋅∂∂+∂∂是可接尺度的变形率张量。

式(14)第2项是为了满足不可压缩的连续方程，当ij S 收缩时(ij S =0)等式两边可以相等。

将亚格子应力的涡粘模型公式(14)代入到(13)式中，变形得
)])([()3(i
j i i t i kk i i i j i x u x u x p x x u u t u ∂∂+∂∂+∂∂++∂∂-=∂∂+∂∂νντρ (16) 0=∂∂i
x u i (17) 4.4.2 Smagorinsky 模型
Smagorinsky 模型是由Smagorinsky 于1963年提出来的，该模型是第一个亚格子模型。

广泛用于大涡模拟中的涡粘模型认为亚格子应力的表达式如下：
(18) 式中(1/2)[(/)(/)]ij i j j i S u x u x =⋅∂∂+∂∂是可接尺度的变形率张量，T ν是涡粘系数。

1963年Smagorinsky 定义了涡粘系数：
ij T kk ij ij S ντδτ231-=-
(19)
式中2/1)
2(ij ij S S S =是变形率张量的大小，∆是过滤尺度，C S 无量纲参数，称为
Smagorinsky 系数。

4.4.3 动态亚格子模式
1991 年, []
8Germano 提出了动态亚格子模式，该模式以Smagorinsky 模式为基本模型，但克服了Smagorinsky 模式的部分缺陷。

动力模型实际上是动态确定亚格子涡粘模型的系数。

动力模型需要对湍流场做两次过滤,一次是细过滤,细过滤后再做一次粗过滤。

通过在网格尺度和检验滤波器尺度条件下计算得到的应力差来确定应力模型系数,使模型系数成为空间和时间的函数,从而避免了在模拟过程中对系数进行调节。

因此比Smagorinsky 模式所采用的固定系数值更加合理。

4.4.4 相似性模式
1980 年Bardina 提出了尺度相似模式。

该模式假定从大尺度脉动到小尺度脉动的动量输运主要由大尺度脉动中的最小尺度脉动来产生,并且过滤后的最小尺度脉度速度和过滤掉的小尺度脉动速度相似。

通过二次过滤和相似性假定可以导出亚格子应力表达式。

采用这种模式能正确预测墙壁面附近的渐近特性,但预测各向不均匀的室内空气复杂流动准确性较差。

4.4.5 混合模式
混合模式是将尺度相似模式和Smagorinsky 模式叠加来确定亚格子应力。

这种模式既有和实际亚格子应力良好的相关性,又有足够的湍动能耗散。

4. 5 优点
(1)能够描述小尺度湍流流动，但是计算量远小于DNS ，在科学研究和工程应用上都显示出良好的发展前景。

(2)用非均匀网格能够使网格数达到最少，节省计算资源，同时又能够保证足够的计算精度。

(3)网格尺度比湍流尺度大，可以模拟湍流发展过程的一些细节。

S C S T 2)(∆=ν
(4)相较于RANS 方法，LES 可以模拟更多的湍流大尺度运动，LES 所用的湍流亚网格应力模型受边界的几何形状和流动类别的影响小，比 RANS 方法所用的Reynolds 应力更具普适性。

4.6 缺点
(1)小涡模型网格节点的划分极密集,需要庞大的计算机存储能力;
(2) 大量数据处理和非线性偏微分方程的求解需要高速数值处理能力;
(3)仅用于比较简单的剪切流运动及管流。

(4) 由于实际湍流极其复杂, 数值模拟仍需要非常可观的计算时间和实验经费。

5. LES,DNS,RANS 三种模拟模型计算量比较
LES, DNS,RANS 三种模拟模型中DNS 的计算量最大，LES 的计算量介于另外两者之间，而RANS 的计算量最小。

影响计算量的因素有三个：网格数量、流场的时间积分长度（与计算时间长度有关）和最小旋涡的时间积分长度（与时间步长有关），其中网格数量是重要因素。

直接数值模拟（DNS ）中为了得到湍流问题足够精确的解，要求能够数值求解所有旋涡的运动，因此要求网格的尺度和最小旋涡的尺度相当，即使采用子域技术，其网格规模也是巨大的。

为了求解各个尺度旋涡的运动，要求每个方向上网格节点的数量与34Re 成比例，考虑一个三维问题，网格节点的数量与94Re 成比例。

一般的估计如下:湍流中包含许多尺度不同的涡，为能模拟最小涡的运动，计算网格的分辨率应足以分辨最小尺度的涡，后者以Kolmogorov 定义的内尺度3=v ηε（）为代表。

而计算区域的尺寸应足以容纳最大尺度的涡，最大涡的尺度为L 。

因此在一个空间方向上的网点数目至少应与L η同量阶，而根据统计理论知道这个比值
3432~L L
R R λη或 （20）
于是整个三维空间所需的网点总数至少为
39942~~L L N R R λη⎛⎫ ⎪⎝⎭或 （21）
此数字也正是按非线性动力系统理论所估计的湍流的吸引子维数的上确界。

计算所需的内存容量应与此数成正比。

另一方面计算的时间步长应小于最小涡的时间尺度'η，而总的计算时间应大于最大涡的特征时间'L u ，因此需要计算的步数应不少于L η.如假设每一时间步长的计算工作量，即使按最低限估计，
与N 成正比，则总的计算工作量至少也要正比于3L R 或6R λ。

假如对每一时间步的每
一网点需执行100条机器指令，则对一个510L R =的湍流问题，就需执行总共约1710条指令。

这意味着在一个计算速度为每秒一亿次的超级计算机上也要运行约30年。

如此巨大的计算工作量即使对当今世界上最大的计算机也是不可接受的。

据Kim ,Moin &Moser 研[]3
究,即使模拟Re 仅为3300 的槽流,所用的网点数N 就约达到了6210⨯ ,在向量计算机上进行了250 h 。

在现有的计算机能力的限制下，即使在少数拥有世界最大的超级计算机的科学大国，目前也只能计算中等以下雷诺数且有简单几何边界的湍流流[]2
动。

湍流大涡数值模拟与直接数值模拟相比节省很大的计算量。

湍流大涡数值模拟将湍流的大尺度脉动和小尺度脉动分开,对大尺度结构进行直接数值模拟,通过建立亚格子尺度(亚格子尺度) 模型来模拟小尺度脉动的作用。

理想的湍流直接数值模拟需要包含所有尺度的湍流脉动,一般最小的脉动尺度等于Kolmogorov 耗散尺度η,流动的最大尺度L 由流动的几何条件确定。

直接数值模拟的一维网格数应为: DNS ~N L η,而大涡数值模拟的一维网格数为: ~LES N L ∆ 可以节省
网格数()()()()33331DNS LES DNS N N N η⎡⎤-=-∆⎣⎦
,如果过滤尺度等于2 倍柯氏耗散尺度的话,就可以比DNS 节省87. 5 %的网格。

这里我们可以看到完全的湍流直接数值模拟中,绝大部分的计算量花费在耗散尺度中,对于高雷诺数流动,这是很不经济的计[]4
算。

雷诺时均方程法先将紊流中的物理量如速度、浓度等分成扰动量及平均量,再利用对控制方程作时间平均,同时采用紊流模型仿真紊流的效应，因此大大降低了计算量,但其结果受紊流模型的影响很大。

参考文献：
[1] 张兆顺. 湍流[M] . 北京:国防工业出版社,2002
[2] 是勋刚. 湍流直接数值模拟进展与前景.北京:北京大学力学系,1992
[3] 路明,孙西欢,李彦军,范志高湍流数值模拟方法及其特点分析.河北建筑科技学院学报,2006
[4] 崔桂香,许春晓,张兆顺. 湍流大涡数值模拟进展[J ] .空气动力学学报,2004 ,22 (2) :121 - 129
[5] 成水燕,李杰,李少飞基于两方程湍流模型的N-S方程数值计算研究.弹箭与制导学报2006,26(1)
[6] 王玲玲. 大涡模拟理论及其应用综述[J]. 河海大学学报, 2004, 32(3):261~265
[7] 肖红林，罗纪生. 大涡模拟中亚格子模型的该进及其在槽道湍流中的应用[J]. 航空动力学报，2007，
22(4):583~587
[8] GERMANO M ,POMELLI U ,et al. A dynamic subgrid -scale eddy viscosity model [J ] .
Phy Fluids ,1991 , (3) :1760- 1765.。