大气边界层在高雷诺数下的模拟结果

大气边界层在高雷诺数下的模拟结果

3.1 不同亚格子模型对大气边界层数值模拟结果的影响

在不可压缩湍流中，亚格子涡粘模型采用分子粘性的形式[19]，即：

123ij ij t ij kk

S τνδτ-=+

(1) 在公式(1)中，t ν为亚格子涡粘系数， _ij S 是可解尺度的变形率张量。可以看出，在涡粘模型中，加上涡粘系数后，就可以使用纳维斯托克斯方程直接数值计算，其中，涡粘系数是需要对其封闭的，可以有不同的形式。

将亚格子涡粘模型(1)代入到大涡数值模拟的控制方程中，可以得到大涡数值模拟控制方程[12]：

_()()3i j j i kk j i j j i u u u u u p t x x x x x τνν---

--∂∂∂∂∂∂+=-+++∂∂∂∂∂∂ ， 1Re t νν=+ (2) _

0i i

u x ∂=∂ (3) 可以看出，代入亚格子涡粘模型和涡扩散模型后，湍流控制方程的变化只是在分子扩散系数上加亚格子涡粘系数，下面给出亚格子涡粘系数的两种形式。

3.1.1 Smagorinsky 模式

在大涡模拟中将雷诺平均混合长度模式推广，得到原始的Smagorinsky 模式，混合长度模式的涡粘公式[20]为：

'2

t u u l l y ν∂<>∝∝∂ 令l =∆，那么在二维空间中有：

2()2ij ij u S S y

∂<>=<><>∂ 将上式推广到三维空间中可得到如下关系： 1

2

2(2)t ij ij l S S ν∝<><>

同样，令l =∆，将平均运算改为过滤运算，并引入模型常数后，可将亚格子涡粘系数写成以下形式：

__222t S ij ij C S S ν=∆<> (4)

在实际模拟计算中发现模型常数明显偏大,Smagorinsky 模式的致命缺点是耗散比较大，下面利用近壁阻尼函数，对Smagorinsky 系数进行修正，从而减少耗散[19]。

对Smagorinsky 系数S C 做如下修正：

0.18[1exp(/)]s C y A ++=--

式中/y yu τν+=，26A +=，[1exp(/)]f y A ++=--称为近壁阻尼函数。

3.1.2 动力模式

涡粘[13]和涡扩散模式在工程计算中有普遍的应用。但因其导出过程以均匀各向同性湍流为基础，实际流动中的复杂湍流并不满足这些条件，所以通过对亚格子模型的改进，使这种模式适应局部湍流的结构，下面介绍的动力模式是经过改进的亚格子涡粘模型。

动力模式在模拟过程中使用了两次过滤，假定过滤过程是线性的，且有12∆>∆，由1∆、2∆过滤产生的亚格子应力和连续两次过滤产生的亚格子应力

可得到如下关系：

13ij kk ij D ij L L C M δ-=

因为上式为超定方程，所以不能直接计算系数D C ,但有以下几种方法来解决超定

问题[14]。第一个为变形率张量收缩法，即将上式两边同乘以可解尺度的变形率张量。但是实际计算表明，用这种方法计算的模式系数很不规则，而且计算过程中的稳定性较差；第二种为最小误差法，令这个超定方程两边的平方差最小，即：

2

103ij kk D ij D L L C M C ⎧⎫∂--=⎨⎬∂⎩⎭ 可得到：

ij ij D ij ij M L C M M = (5)

在式(5)中，ij L 和流动的形态有关，ij M 和模式的形式有关。最小误差法比变形

率张量收缩法已经有很大改进，但是仍有缺陷：①模式系数D C 可能会出现负值，

导致计算结果发散；②分母可能很小，也导致计算结果发散。为了克服计算上的困难，采用平均系数法，对(5)式右端的分子和分母分别求系综平均，得到模式系数为：

ij ij D ij ij M L C M M <>=<> (6)

3.1.3 大涡模拟的验证

将大涡模拟应用到低雷诺数的三维周期性槽道流中，使用的亚格子模型为动力模式，雷诺数的值为4400，因为大涡模拟不直接计算近壁区域，所以增大了网格的分辨率[15]，模拟过程中使用的网格点数为：365，计算域和直接数值模拟相同:02,02,0X Y Z ≤≤∏≤≤≤≤∏，0.01t ∆=。下图为得到的壁面应力随时间的变化情况：

图3.1 大涡模拟壁面应力图

图3.1为用大涡模拟计算的壁面应力随时间的演化过程。可见，在无量纲时间200之后，数据趋于稳定，所以在计算平均速度以及脉动和雷诺应力时，从20000步开始取数据，接下来画出流向速度随壁面高度的变化，处理方法同直接数值模拟：

图3.2 大涡模拟流向速度图

图3.2为用大涡模拟计算的流向速度随时间的演化过程。因为未加壁函数，所以第一个网格点在近壁区域，大涡模拟取得网格点比直接数值模拟少，得到的结果没有直接数值模拟精确，计算结果与正确的结果之间出现了一定的偏差，但是，因为大涡模拟所取的网格点比较稀疏，加之壁面对结果的影响比较大，可看出，得到的结果是在误差范围内的。下图为脉动''u u <>、''v v <>、''w w <>

和雷诺应力''u v <>随壁面高度的变化，在平均上与速度的处理方式一样，先对

同一个时间点同一高度上的值进行平均，然后在各个时间上平均，实线代表模拟计算得到的结果：

图3.3 大涡模拟流向脉动图 图3.4 大涡模拟法向脉动图

图3.5 大涡模拟展向脉动图 图3.6 大涡模拟雷诺应力图 图3.3、图3.4、图3.5、图3.6分别为用大涡模拟计算的流向、法向、展向脉动和雷诺应力随时间的演化过程。所以，在大涡模拟中，速度脉动和雷诺应力在误差范围内也得到了相对较好的验证。同时也证明了现有程序对大涡模拟的有效性。

3.1.4 两种亚格子模型的模拟结果

在3.1.1和3.1.2中分别介绍了Smagorinsky 模式和动力模式，并且利用近壁阻尼函数对Smagorinsky 模式进行了改进，动力模式采用涡粘模式表达式，改进之后的模型在减少计算量的基础上有效的模拟了壁面对湍流边界层的影响，在模拟湍流边界层的过程中雷诺数增加到了15000，网格数为365，计算域不变 ，0.005t ∆=，缩短了时间间隔，首先给出两种亚格子模型下流向速度随壁面高度的变化趋势，并将结果进行对比：

图3.7 不同亚格子模型模拟流向速度图