文章快速检索 高级检索
0

页岩气专题论文

引用本文 [复制中英文]

秦泽聪, 方乐. 一种改进的均匀各向同性湍流初始化方法1)[J]. 力学学报, 2016, 48(6): 1319-1325. DOI: 10.6052/0459-1879-16-180.
[复制中文]
Qin Zecong, Fang Le. AN IMPROVED METHOD FOR INITIALIZING HOMOGENEOUS ISOTROPIC TURBULENT FLOWS1)[J]. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(6): 1319-1325. DOI: 10.6052/0459-1879-16-180.
[复制英文]

基金项目

1) 国家自然科学基金资助项目(11202013,11302012,51420105008)

通讯作者

2) 方乐, 副教授, 主要研究方向:流体力学.E-mail:le.fang@zoho.com

文章历史

2016-07-01 收稿
2016-09-04 录用
2016-09-21网络版发表
一种改进的均匀各向同性湍流初始化方法1)
秦泽聪*,†, 方乐*2)     
*. 北京航空航天大学中法工程师学院, 北京 100191
†. 法国里昂中央理工大学流体与声学实验室, 里昂 69130, 法国
摘要: 均匀各向同性湍流是一种最简单的湍流理想状态,也是湍流基础理论研究的最重要对象之一.为了用数值方法产生均匀各向同性湍流场,一般采用Rogallo提出的方法在谱空间生成初始场,然后再转换到物理空间.研究表明,由该方法生成的初始湍流场在3个棱向上呈各向异性,在结构函数和速度概率密度分布上均有体现.尽管在初始场样本很多时,这种各向异性可以在平均意义上消除,但作为数值模拟采用的单个流场则波动较大,不利于在实际计算中作为单个初始场生成各向同性湍流.在此基础上提出一种改进的Rogallo初始化方法,称为模量平均法,将Rogallo方法在3个轴向分别进行,并进行模量平均,最后采用能谱进行模量控制.这种方法可以一方面保持初始场能谱,另一方面减小单个流场的各向异性波动,以产生各向同性程度更佳的单个初始场.在统计意义上,新方法可以分别将结构函数和速度概率密度的相对标准差减小约10%.
关键词: 均匀各向同性湍流    Rogallo初始化方法    结构函数    概率密度分布    
引言

均匀各向同性湍流是一种理想的湍流简化形态[1-2].为了在数值上生成均匀各向同性湍流,目前研究者广泛采用的做法是在三维空间内引入周期性边界条件[3-10],将无限大空间近似为无穷多个周期性的立方体.这种做法其实是由于数值分辨率的限制,必须在谱空间进行尺度截断而引入的.小尺度的截断与数值误差、能量积累等问题有关[11-12],而周期性条件则是一种最常见的大尺度的截断.这种大尺度截断会引起一系列的问题,包括间歇性[13]、各向异性等问题.显然,如果由于引入三维周期性条件带来额外的各向异性,势必会影响我们希望生成的“各向同性”湍流场的质量.

文献[14]对此问题已经进行了初步的讨论,将湍流的各向异性问题归结为张量各向异性和标量方向各向异性两种.前者指的是二阶纵向结构函数Dll (下标l表示结构函数纵向分量)和横向结构函数Dnn (下标n表示结构函数横向分量)在不同矢径长度r下需满足各向同性关系式[15-18]

$ {D_{{\rm{nn}}}}(r) = {D_{{\rm{ll}}}}(r) + \frac{r}{2}\frac{{{\rm{d}}{D_{{\rm{ll}}}}(r)}}{{{\rm{d}}r}} $ (1)

这在一些数值模型验证工作中曾被作为各向同性检测方法使用[19-20].后者则指对于标量两点统计量,其统计值应与方向无关,只与两点距离有关.前期工作已经证明[14],采用三维周期性条件必定会在大尺度引起棱与对角线之间的标量方向各向异性,其引起的各向异性差异可达10%.这种各向异性可以通过引入六边形谱方法来尽量减小[21].

本文将关注采用三维周期性条件时初始场3个棱方向的各向异性问题.此各向异性问题不仅在二阶纵向结构函数中体现,而且甚至可在各个速度分量的概率密度分布中体现,这是一种由于样本波动较大而引起的统计问题.本文提出一种改进方法,以减小样本波动,生成各向同性程度更好的初始湍流场.

1 Rogallo初始场生成方法及其棱向各向异性 1.1 Rogallo初始场生成方法

在进行计算机数值模拟时,均匀各向同性湍流场的流场区域一般被设计为一个棱长为2π的立方体形状.这种流场普遍使用Rogallo[22]在1981年提出的方法进行初始化.该方法的核心思想是在谱空间中构造一个同时满足连续性方程以及给定能谱的流场

$ {k_1}\hat U(\mathit{\boldsymbol{k}}) + {k_2}\hat V(\mathit{\boldsymbol{k}}) + {k_3}\hat W(\mathit{\boldsymbol{k}}) = 0 $ (2)
$ \mathop{{\int\!\!\!\!\!\int}\mkern-21mu \bigcirc}\limits_{A(k)} {\mathit{\boldsymbol{\widehat u}}\left( \mathit{\boldsymbol{k}} \right) \cdot {{\mathit{\boldsymbol{\widehat u}}}^*}\left( \mathit{\boldsymbol{k}} \right)} {\rm{d}}A(k) = E(k) $ (3)

其中,ki$ \hat {U} $, $ \hat {V} $, $ \hat {W} $分别为谱空间波数矢量k和速度矢量$ \boldsymbol{\widehat u} $在3个棱方向上的分量,$ \boldsymbol{\widehat u} $*为速度矢量$ \boldsymbol{\widehat u} $的共轭矢量,A(k)为谱空间中半径为k(波数矢量k的模长)的球面积分域,E(k)为给定能谱.连续性方程要求谱空间中速度矢量与波数矢量垂直[23].因此,Rogallo方法根据每一个波数矢量k建立局部坐标基(e'1, e'2, e'3),并令e'3矢量与k矢量平行,即

$ {{\mathit{\boldsymbol{e'}}}_3}(\mathit{\boldsymbol{k}}) = \frac{{{k_1}}}{k}{\mathit{\boldsymbol{e}}_1} + \frac{{{k_2}}}{k}{\mathit{\boldsymbol{e}}_2} + \frac{{{k_3}}}{k}{\mathit{\boldsymbol{e}}_3} $ (4)

其中,e1, e2, e3分别为坐标基矢量.则速度矢量$ \boldsymbol{\widehat u} $必须处于由e'1, e'2张成的平面上,写为

$ \mathit{\boldsymbol{\hat u}}(\mathit{\boldsymbol{k}}) = \alpha (k){{\mathit{\boldsymbol{e'}}}_1}(\mathit{\boldsymbol{k}}) + \beta (k){{\mathit{\boldsymbol{e'}}}_2}(\mathit{\boldsymbol{k}})\; $ (5)

式中,α(k)与β(k)为速度矢量在方向e'1e'2方向上的分量.由于速度场还需要满足能谱方程,因此,两个分量的最终形式为

$ \alpha (k) = \sqrt {\frac{{E\left( k \right)}}{{4{\rm{\pi }}{k^2}}}} \exp \left( {{\rm{i}}{\theta _1}} \right)\cos \phi $ (6)
$ \beta (k) = \sqrt {\frac{{E\left( k \right)}}{{4{\rm{\pi }}{k^2}}}} \exp \left( {{\rm{i}}{\theta _2}} \right)\sin \phi $ (7)

式中,θ1, θ2φ都是在[0, 2π)区间上均匀分布的随机标量.令

$ {{\mathit{\boldsymbol{e'}}}_1}(\mathit{\boldsymbol{k}}) = \frac{{{{\mathit{\boldsymbol{e'}}}_3}(\mathit{\boldsymbol{k}}) \times {\mathit{\boldsymbol{e}}_3}(\mathit{\boldsymbol{k}})}}{{\left| {{{\mathit{\boldsymbol{e'}}}_3}(\mathit{\boldsymbol{k}}) \times {\mathit{\boldsymbol{e}}_3}(\mathit{\boldsymbol{k}})} \right|}} $ (8)
$ {{\mathit{\boldsymbol{e'}}}_2}(\mathit{\boldsymbol{k}}) = \frac{{{{\mathit{\boldsymbol{e'}}}_3}(\mathit{\boldsymbol{k}}) \times {{\mathit{\boldsymbol{e'}}}_1}(\mathit{\boldsymbol{k}})}}{{\left| {{{\mathit{\boldsymbol{e'}}}_3}(\mathit{\boldsymbol{k}}) \times {{\mathit{\boldsymbol{e'}}}_1}(\mathit{\boldsymbol{k}})} \right|}} $ (9)

得到谱空间中的Rogallo初始场

$ \mathit{\boldsymbol{\hat u}}(\mathit{\boldsymbol{k}}) = \frac{1}{{k\sqrt {k_1^2 + k_2^2} }}\left( {\begin{array}{*{20}{c}} {\alpha (k)k{k_2} + \beta (k){k_1}{k_3}}\\ {\beta (k){k_2}{k_3}-\alpha (k)k{k_1}}\\ {-\beta (k)(k_1^2 + k_2^2)} \end{array}} \right) $ (10)
1.2 Rogallo初始场的棱向各向异性

采用一个典型的算例来说明采用Rogallo方法生成的初始场的棱向各向异性问题. 图 1所示为某个空间网格点数为2563的Rogallo初始场的二阶纵向结构函数Dll在3个棱方向上的曲线(两点距离Y用网格尺度进行无量纲化;结构函数用脉动速度均方根的平方进行无量纲化).可以观察到,3条曲线在大尺度处差异较大,显示该初始场在棱之间具有明显的各向异性.需要说明的是,虽然在这里只给出了某个初始场的统计结果,但这种现象不是个例,在实际计算中非常普遍.

图 1 单个Rogallo初始场的各棱向二阶纵向结构函数Dll Figure 1 Second-order longitudinal structure function Dll of a Rogallo initial field on axis directions

为了对比,我们也同时计算了采用Rogallo方法生成的初始场各棱向速度分量U, V, W的概率密度函数(possibility density function,PDF),如图 2所示(横坐标速度用脉动速度的均方根进行无量纲化).同样可以看到,3个棱向速度分量的PDF之间存在明显的差异.由于PDF和一维能谱之间存在对应关系, 这种差异显示了3个棱向湍动能之间分配的各向异性,也正是我们希望在具体应用中尽量减小的.下面我们从统计意义上给这种差异以定量的描述.

图 2 单个Rogallo初始场的各棱向速度分量概率密度 Figure 2 Axis direction velocity component possibility density of a Rogallo initial field
2 针对Rogallo方法统计性质的改进方案 2.1 Rogallo初始场的统计性质

从结构函数的各向异性结果中可以看出,棱向各向异性明显的区域位于大尺度区.这可能是由于在单个初始场中大尺度的流动结构数量太少而引起的样本缺乏所导致的.为了证实这一结论,需要增加样本数以研究初始场之间的统计性质.通过改变式(6)和式(7)中随机序列θ1, θ2φ的随机种子,可以得到满足相同能谱的不同初始场.通过计算机程序生成1 024个基于Rogallo方法的、空间网格点数均为2563的初始场样本,每个初始场的棱向Dll(r)以及速度3个棱向分量的PDF均存在类似图 1图 2所示的差异.但如果对这1 024个Dll标量场及PDF曲线做求和平均处理,则可以观察到这些差异几乎被完全消除,结果如图 3图 4所示.这说明Rogallo初始化方法本身并不会在平均意义上引起棱向各向异性.尽管如此,实际计算中使用的并不是上千个初始场的平均,而只是单个的初始场,仍然不希望每个初始场的波动太大.因此,有必要对初始场之间的波动进行统计计算.需要强调的是,在这里一个“样本”指的是一个完整的初始场,而并非一个离散点,这是由于实际计算中我们总是采用单个初始场,因而对这种单个样本的波动特性具有较高的要求.总之,在这里将单个初始场视为样本是较合理的做法.

图 3 多个Rogallo初始场平均后的各棱向二阶纵向结构函数Dll(无量纲方法同图 1) Figure 3 Averaged second-order longitudinal structure function Dll of Rogallo initial fields on aixs directions (the non-dimensialization information can be found in Fig. 1)
图 4 多个Rogallo初始场平均后的各棱向速度分量概率密度(无量纲方法同图 2) Figure 4 Averaged axis direction velocity component possibility density of Rogallo initial fields (the non-dimensialization information can be found in Fig. 2)

为了定量地表征这些初始场样本的波动大小,定义相对标准差

$ {\sigma _r}(a(t)) = \sqrt {\frac{{\frac{1}{N}\sum\limits_{n = 1}^N {{{({a_n}(t)-\mu )}^2}} }}{{{\mu ^2}}}} $ (11)

式中,N=1 024为初始场的总样本数,an(t)表示所统计的对象a(t)的第n个样本,t为相应的自变量,μa(t)的平均值,即

$ \mu = \frac{1}{N}\sum\limits_{n = 1}^N {{a_n}(t)} $ (12)

如果我们讨论的是用球坐标表示的二阶纵向结构函数Dll(r, θ, φ)在3个棱向上的曲线各自的样本波动情况,则a(t)分别为Dll, θ=π/2, φ=0(r), Dll, θ=π/2, φ=π/2(r), Dll, θ=0, φ=0(r);如果讨论的是3个棱向速度分量PDF各自的样本波动情况,则a(t)分别为U, V, W的PDF.这6个对象的相对标准差统计结果分别如图 5图 6所示.从图 5可以看到,Dll(r)的样本波动程度在大尺度更明显,达到约5%. 图 6表明当速度分量的值越靠近平均值0时,其相对波动程度越弱;而在脉动速度较大的区域,由于速度在此处概率密度很低,样本之间的差异非常大,其相对标准差的统计值也因此很大,可达2 000%左右.这种异常大的波动可能会导致额外的间歇性,与Brun等[13]的工作可能有一定联系,此处不做深入讨论.本文中主要关注湍流能量的主要存在区域,即均值附近的区域,相应的可画出此区域U, V, W的PDF相对标准差,如图 7所示.可以观察到此区域的相对标准差大约为3%.

图 5 多个Rogallo初始场各棱向二阶纵向结构函数的相对标准差(两点距离用网格尺度进行无量纲化) Figure 5 Relative standard deviations of the second-order longitudinal structure function of Rogallo initial fields on axis directions (the two-point distance is non-dimensionalized by the grid length )
图 6 多个Rogallo初始场各棱向速度分量概率密度的相对标准差(无量纲方法同图 2) Figure 6 Relative standard deviations of axis direction velocity component possibility density (PDF) of Rogallo initial fields (the non-dimensialization information can be found in Fig. 2)
图 7 多个Rogallo初始场各棱向速度分量概率密度的相对标准差(无量纲脉动速度绝对值小于0.4部分) (无量纲方法同图 2) Figure 7 Relative standard deviations of axis direction velocity component possibility density (PDF) of Rogallo initial fields (only the central zone where the non-dimensionalized velocity absolute value is smaller than 0.4) (The non-dimensialization information can be found in Fig. 2)
2.2 模量平均法

为了尽量减小初始场样本之间的波动,得到各向同性程度更好的初始场,需要对Rogallo方法进行分析.事实上,式(10)表明谱空间速度矢量$ \mathit{\boldsymbol{\hat u}}(\mathit{\boldsymbol{k}}) $的第1、第2个分量与α(k)和β(k)均显式相关,但第3个分量只与β(k)显式相关,这导致$ \mathit{\boldsymbol{\hat u}}(\mathit{\boldsymbol{k}}) $的3个分量在谱空间是形式上统计不对称的.经过傅里叶逆变换之后,物理空间的棱向速度分量U, V, W之间也因此存在不对称性.

为了使$ \mathit{\boldsymbol{\hat u}}(\mathit{\boldsymbol{k}}) $三个分量之间在谱空间统计对称,可以对每个k矢量分别令ke'3, e'1, e'2平行,同时生成3个随机速度场

$ {{\mathit{\boldsymbol{\hat u}}}_1}(\mathit{\boldsymbol{k}}) = \frac{1}{{k\sqrt {k_1^2 + k_2^2} }}\left( {\begin{array}{*{20}{c}} {{\alpha _1}(k)k{k_2} + {\beta _1}(k){k_1}{k_3}}\\ {{\beta _1}(k){k_2}{k_3}-{\alpha _1}(k)k{k_1}}\\ {-{\beta _1}(k)(k_1^2 + k_2^2)} \end{array}} \right) $ (13)
$ {{\mathit{\boldsymbol{\hat u}}}_2}(\mathit{\boldsymbol{k}}) = \frac{1}{{k\sqrt {k_2^2 + k_3^2} }}\left( {\begin{array}{*{20}{c}} {-{\beta _2}(k)(k_2^2 + k_3^2)}\\ {{\alpha _2}(k)k{k_3} + {\beta _2}(k){k_2}{k_1}}\\ {{\beta _2}(k){k_3}{k_1}-{\alpha _2}(k)k{k_2}} \end{array}} \right) $ (14)
$ {{\mathit{\boldsymbol{\hat u}}}_3}(\mathit{\boldsymbol{k}})\frac{1}{{k\sqrt {k_1^2 + k_2^2} }}\left( {\begin{array}{*{20}{c}} {{\beta _3}(k){k_1}{k_2}-{\alpha _3}(k)k{k_3}}\\ {-{\beta _3}(k)(k_3^2 + k_1^2)}\\ {{\alpha _3}(k)k{k_1} + {\beta _3}(k){k_3}{k_2}} \end{array}} \right) $ (15)

式中,αi(k)和βi(k)均是形如式(6)、式(7)的随机数,但各自采用不同的随机序列.接着对这3个矢量场作如下平均,得到改进的湍流初始场

$ \mathit{\boldsymbol{\hat u}}(\mathit{\boldsymbol{k}}) = \sqrt {\frac{{E(k)}}{{4\pi {k^2}}}} \frac{{({{\mathit{\boldsymbol{\hat u}}}_1}(\mathit{\boldsymbol{k}}) + {{\mathit{\boldsymbol{\hat u}}}_2}(\mathit{\boldsymbol{k}}) + {{\mathit{\boldsymbol{\hat u}}}_3}(\mathit{\boldsymbol{k}}))}}{{|{{\mathit{\boldsymbol{\hat u}}}_1}(\mathit{\boldsymbol{k}}) + {{\mathit{\boldsymbol{\hat u}}}_2}(\mathit{\boldsymbol{k}}) + {{\mathit{\boldsymbol{\hat u}}}_3}(\mathit{\boldsymbol{k}})|}} $ (16)

我们称这种改进方法为模量平均法.容易看出式(16)能同时满足谱空间中的连续方程(2)和能谱方程(3),且在谱空间3个方向上统计对称.如果把式(16)代入能谱方程(3)的左边部分,给定能谱可自动满足

$ \begin{array}{*{20}{c}} {\mathop{{\int\!\!\!\!\!\int}\mkern-21mu \bigcirc}\limits_{A(k)} {\mathit{\boldsymbol{\widehat u}}\left( \mathit{\boldsymbol{k}} \right) \cdot {{\mathit{\boldsymbol{\widehat u}}}^*}\left( \mathit{\boldsymbol{k}} \right)} {\rm{d}}A(k) = \mathop{{\int\!\!\!\!\!\int}\mkern-21mu \bigcirc}\limits_{A(k)} {\frac{{E(k)}}{{4{\rm{ \mathsf{ π} }}{k^2}}}{\rm{d}}A(k)} = }\\ {\frac{{E(k)}}{{4{\rm{ \mathsf{ π} }}{k^2}}}\mathop{{\int\!\!\!\!\!\int}\mkern-21mu \bigcirc}\limits_{A(k)} {{\rm{d}}A(k)} = E(k)} \end{array} $ (17)

图 8是实际计算结果,可以看到式(16)生成的初始场能满足给定能谱.

图 8 基于模量平均法生成的初始场能谱与给定能谱的对比图 Figure 8 Comparison between the given spectrum and the spectrum calculated from an initial field based on the modulus-averaging method

为了观察该方法的统计性质,采用类似1.3节的做法,我们通过修改随机序列的种子生成1 024个基于模量平均法的初始场,并相应地计算3个棱向Dll(r)的相对标准差以及U, V, W的PDF相对标准差,结果分别如图 9图 10 (无量纲脉动速度绝对值小于0.4的部分)所示.通过图 5图 9图 7图 10的对比,可以看出各物理量的样本波动程度都有所减弱.定量的比较如图 11图 12 (无量纲速度脉动值小于0.4的部分)所示.可以看到,在大多数区域,采用模量平均法能将初始场的相对标准差降低约10%.

图 9 Comparison between the given spectrum and the spectrum calculated from an initial field based on the modulus-averaging method Figure 9 Relative standard deviationsof the second-order longitudinal structure function on axes directions of initial fields generated by the modulus-averaging method (The non-dimensialization information can be found in Fig. 5)
图 10 基于模量平均法生成的多个初始场各棱向速度分量概率密度的相对标准差(无量纲方法同图 2) Figure 10 Relative standard deviations of axis direction velocity component possibility density (PDF) of initial fields generated by the modulus-averaging method (The non-dimensialization information can be found in Fig. 2)
图 11 多个初始场各棱向二阶纵向结构函数的相对标准差,采用不同初始化方法时的相对差异(定义为(模量平均法Rogallo方法)-(Rogallo方法))(无量纲方法同图 5) Figure 11 The relative difference between relative standard deviations (RSD) of second-order longitudinal structure function on axis directions by different initialization methods (The definition is (modulus-averaging method-Rogallo method)-(Rogallo method)) (The non-dimensialization information can be found in Fig. 5)
图 12 多个初始场各棱向速度分量概率密度的相对标准差,采用不同初始化方法时的相对差异(定义为(模量平均法Rogallo方法)-(Rogallo方法))(无量纲方法同图 2) Figure 12 The relative difference between relative standard deviations (RSD) of axis direction velocity component possibility density (PDF) by different initialization methods (The definition is (modulus-averaging method-Rogallo method)-(Rogallo method)) (The non-dimensialization information can be found in Fig. 2)
3 结论

为了生成均匀各向同性湍流场,我们要求初始场应尽可能接近真实的均匀各向同性湍流.尽管某些表征湍流统计状态的统计量(如表征湍流非平衡性的速度梯度扭率[24-26]和四阶矩[27],表征湍流时空特性的时空相关[28-29])很难在初始场中加以构造,但至少对于结构函数、速度概率密度这种与能谱相关的简单的统计量,初始场应当尽量各向同性,这应当是评价初始场质量的一个基本要求.初始场将能显著影响湍流非均衡发展的时间[30].

目前被广泛采用的初始场生成方法由Rogallo提出.在实际应用中,Rogallo方法生成的初始场并不能很好保证结构函数和速度概率密度的各向同性,如果将每个初始场视为一个统计样本,则样本波动较大.在Rogallo方法的基础上我们提出一种模量平均法.与Rogallo方法相比,新方法能分别将结构函数和速度概率密度的相对标准差减小约10%,在一定程度上有所改进.因此我们建议,在数值计算中可以采用模量平均法来取代Rogallo方法生成初始场.


参考文献
1 张兆顺, 崔桂香, 许春晓. 湍流大涡数值模拟的理论与应用 . 北京: 清华大学出版社, 2008. ( Zhang Zhaoshun, Cui Guixiang, Xu Chunxiao. Theory and Application of Large Eddy Numerical Simulation . Beijing: Tsinghua University Press, 2008. (in Chinese) )
2 Pope SB. Turbulent Flows . New York: Cambridge University Press, 2000.
3 李瑞霞, 柳朝晖, 贺铸, 等. 各向同性湍流内颗粒碰撞率的直接模拟研究. 力学学报, 2006, 38 (1) : 25-32. ( Li Ruixia, Liu Zhaohui, He Zhu, et al. Direct numerical simulation of inertial particle collisions in isotropic turbulence. Chinese Journal of Theoretical and Applied Mechanics, 2006, 38 (1) : 25-32. (in Chinese) )
4 李新亮, 傅德薰, 马延文. 可压缩均匀各向同性湍流的直接数值模拟. 中国科学(A辑), 2002, 32 (8) : 716-724. ( Li Xinliang, Fu Dexun, Ma Yanwen. Direct numerical simulation of compressible homogeneous isotropic turbulence. Science in China (Series A), 2002, 32 (8) : 716-724. (in Chinese) )
5 张曙光, 雷磊, 何国威. 均匀各向同性湍流的频率波数能量谱. 力学学报, 2003, 35 (3) : 317-320. ( Zhang Shuguang, Lei Lei, He Guowei. Frequency-wavenumber energy spectra in isotropic turbulence. Acta Mechanica Sinica, 2003, 35 (3) : 317-320. (in Chinese) )
6 赵松年, 胡非. 湍流问题:如何看待"均匀各向同性湍流"?. 中国科学:物理学力学天文学, 2015, 45 (2) : 024701. ( Zhao Songnian, Hu Fei. Turbulence question:How do view "the homogeneous and isotropic turbulence"?. Science China Physics, Mechanics & Astronomy, 2015, 45 (2) : 024701. (in Chinese) )
7 方乐, 杨云柯, 王洪涛, 等. 全场离散Tophat过滤操作的快速算法. 数值计算与计算机应用, 2009, 30 (3) : 218-224. ( Fang Le, Yang Yunke, Wang Hongwei, et al. A rapid algorithm for Tophat filter operation in discrete field. Journal on Numerical and Computer Applications, 2009, 30 (3) : 218-224. (in Chinese) )
8 方乐, 崔桂香, 许春晓, 等. 标量湍流的能量传输特性. 计算物理, 2006, 23 (6) : 692-698. ( Fang Le, Cui Guixiang, Xu Chunxiao, et al. Energy transfer in scalar turbulence. Chinese Journal of Computational Physics, 2006, 23 (6) : 692-698. (in Chinese) )
9 Stepanov R, Plunian F, Kessar M, et al. Systematic bias in the calculation of spectral density from a three-dimensional spatial grid. Physical Review E, 2014, 90 (5) : 053309. DOI: 10.1103/PhysRevE.90.053309.
10 Fang L, Cui GX, Xu CX, et al. Multi-scale analysis of energy transfer in scalar turbulence. Chinese Physics Letters, 2005, 22 (11) : 2877. DOI: 10.1088/0256-307X/22/11/042.
11 Cichowlas C, Bonaïti P, Debbasch F, et al. Effective dissipation and turbulence in spectrally truncated Euler flows. Physical Review Letters, 2005, 95 (26) : 264502. DOI: 10.1103/PhysRevLett.95.264502.
12 Bos WJT, Bertoglio JP. Dynamics of spectrally truncated inviscid turbulence. Physics of Fluids, 2006, 18 (7) : 071701. DOI: 10.1063/1.2219766.
13 Brun C, Pumir A. Statistics of Fourier modes in a turbulent flow. Physical Review E, 2001, 63 (5) : 056313. DOI: 10.1103/PhysRevE.63.056313.
14 Qin ZC, Fang L, Fang J. How isotropic are turbulent flows generated by using periodic conditions in a cube?. Physics Letters A, 2016, 380 (13) : 1310-1317. DOI: 10.1016/j.physleta.2016.02.001.
15 Hinze JO. Turbulence . New York: McGraw-Hill International Book Company, 1975.
16 Fang L, Bos WJT, Jin GD. Short-time evolution of Lagrangian velocity gradient correlations in isotropic turbulence. Physics of Fluids, 2015, 27 (12) : 125102. DOI: 10.1063/1.4936140.
17 马威, 方乐, 邵亮, 等. 可解尺度各向同性湍流的标度律. 力学学报, 2011, 43 (2) : 267-276. ( Ma Wei, Fang Le, Shao Liang, et al. Scaling law of resolved-scale isotropic turbulence. Chinse Journal of Theoretical and Applied Mechanics, 2011, 43 (2) : 267-276. (in Chinese) )
18 Fang L, Shao L, Bertoglio JP, et al. An improved velocity increment model based on Kolmogorov equation of filtered velocity. Physics of Fluids, 2009, 21 (6) : 065108. DOI: 10.1063/1.3153911.
19 Gotoh T, Fukayama D, Nakano T. Velocity field statistics in homogeneous steady turbulence obtained using a high-resolution direct numerical simulation. Physics of Fluids, 2002, 14 (3) : 1065-1081. DOI: 10.1063/1.1448296.
20 Cui GX, Zhou HB, Zhang ZS, et al. A new dynamic subgrid eddy viscosity model with application to turbulent channel flow. Physics of Fluids, 2004, 16 (8) : 2835-2842. DOI: 10.1063/1.1762911.
21 乔海军, 李会元. 二维各向同性湍流直接数值模拟的六边形谱方法及GPU实现和优化. 数值计算与计算机应用, 2013, 34 (6) : 147-160. ( Qiao Haijun, Li Huiyuan. Hexagonal spectral method for direct numerical simulation of two-dimensional homogeneous isotropic turbulence and their GPU implementation and optimization. Journal on Numerical Methods and Computer Applications, 2013, 34 (6) : 147-160. (in Chinese) )
22 Rogallo RS. Numerical experiments in homogeneous turbulence. Technical Report 81315, NASA, 1981
23 Lesieur M. Turbulence in Fluids. Springer Science & Business Media, 2012
24 Fang L, Zhu Y, Liu YW, et al. Spectral non-equilibrium property in homogeneous isotropic turbulence and its implication in subgridscale modeling. Physics Letters A, 2015, 379 (38) : 2331-2336. DOI: 10.1016/j.physleta.2015.05.029.
25 Fang L, Bos WJT, Shao L, et al. Time reversibility of Navier-Stokes turbulence and its implication for subgrid scale models. Journal of Turbulence, 2012 (13) : N3.
26 Hearst RJ, Lavoie P. Velocity derivative skewness in fractalgenerated, non-equilibrium grid turbulence. Physics of Fluids, 2015, 27 (7) : 071701. DOI: 10.1063/1.4926356.
27 Fang L, Zhang YJ, Fang J, et al. Relation of the fourth-order statistical invariants of velocity gradient tensor in isotropic turbulence. Physical Review E, 2016, 94 (2) : 023114. DOI: 10.1103/PhysRevE.94.023114.
28 Zhao X, He GW. Space-time correlations of fluctuating velocities in turbulent shear flows. Physical Review E, 2009, 79 (5) : 046316.
29 He GW, Zhang JB. Elliptic model for space-time correlations in turbulent shear flows. Physical Review E, 2006, 73 (5) : 055303. DOI: 10.1103/PhysRevE.73.055303.
30 Vassilicos JC. Dissipation in turbulent flows. Annual Review of Fluid Mechanics, 2015, 47 : 95-114. DOI: 10.1146/annurev-fluid-010814-014637.
AN IMPROVED METHOD FOR INITIALIZING HOMOGENEOUS ISOTROPIC TURBULENT FLOWS1)
Qin Zecong*,†, Fang Le*2)     
*. Sino-French Engineer School, Beihang University, Beijing 100191, China;
†. LMFA, CNRS, Ecole Centrale de Lyon-Université de Lyon, Ecully 69130, France
Abstract: Homogeneous isotropic turbulence (HIT) is one of the simplest ideal turbulence states, and is also one of the most important subjects in basic turbulence theory researches. HIT fields are usually initialized in the spectrum space via the method proposed by Rogallo, and then transformed in physical space. The current paper points out that initial fields thus generated are anisotropic in axis directions of their computational domain, which can be reflected in structure functions and in possibility density distribution of velocity components. Even though such anisotropy will disappear after an average operation of a large number of initial field samples, the anisotropy fluctuation between samples is considerately big, which is not favorable for the establishment of an HIT. Basing on this existing methodology, we then proposed an improved Rogallo method, named the modulus-averaging method, which firstly conducts the Rogallo method in all the 3 axis directions, then carries out a modulus-averaging operation, and finally control the modulus via a given spectrum function. This method can keep the initial filed spectrum and, reduce the anisotropy fluctuation of each single field to generate "more isotropic" initial fields. Statistically, this new method can lower the relative standard deviations of structure functions and velocity possibility density distribution by about 10%.
Key words: homogeneous isotropic turbulence    rogallo initialization method    structure function    possibility density distribution