船舶力学

网格方案	网格总量N（×10⁶）	时间步长（$ \Delta {t}^{+}=\Delta t{U}_{\mathrm{\infty }}/D $）	统计平均间隔（$ \Delta {T}^{+}=\Delta T{U}_{\mathrm{\infty }}/D $）
Grid 1 粗糙网格	2.85	2.5×10⁻³	300
Grid 2 精细网格	3.47	2.1×10⁻³	300
Grid 3 超精细网格	4.12	1.6×10⁻³	300

网格方案	网格总量N（×10⁶）	时间步长（$ \Delta {t}^{+}=\Delta t{U}_{\mathrm{\infty }}/D $）	统计平均间隔（$ \Delta {T}^{+}=\Delta T{U}_{\mathrm{\infty }}/D $）
Grid 1 粗糙网格	2.85	2.5×10⁻³	300
Grid 2 精细网格	3.47	2.1×10⁻³	300
Grid 3 超精细网格	4.12	1.6×10⁻³	300

编号	坐标$ \left(x/D,y/D\right) $	编号	坐标$ \left(x/D,y/D\right) $
1	（0.25, 0.00）	4	（0.25, 0.10）
2	（1.00, 0.00）	5	（1.00, 0.10）
3	（1.75, 0.00）	6	（1.75, 0.10）

编号	坐标$ \left(x/D,y/D\right) $	编号	坐标$ \left(x/D,y/D\right) $
1	（0.25, 0.00）	4	（0.25, 0.10）
2	（1.00, 0.00）	5	（1.00, 0.10）
3	（1.75, 0.00）	6	（1.75, 0.10）

相关研究文献	$ S {t}_{1,{\mathrm{Sub}}} $	$ S {t}_{1} $	$ S {t}_{2} $	$ S {t}_{2,{\mathrm{Mul}}} $
Pereira等^[4]（EXP）			0.55
Pereira等^[5]（EXP）			0.54~0.58	1.24（倍频） 1.72（三倍频）
Chang 等^[8] （LES）	0.185	0.38	0.51
Chang 等^[9] （SST-DES）	/	0.31~0.35	0.51
Rockwell（势流）	/	0.38	/
Ahuja & Mendoza^[22]（理论）	/	0.35	0.505
Howe^[23]（理论）		0.39

相关研究文献	$ S {t}_{1,{\mathrm{Sub}}} $	$ S {t}_{1} $	$ S {t}_{2} $	$ S {t}_{2,{\mathrm{Mul}}} $
Pereira等^[4]（EXP）			0.55
Pereira等^[5]（EXP）			0.54~0.58	1.24（倍频） 1.72（三倍频）
Chang 等^[8] （LES）	0.185	0.38	0.51
Chang 等^[9] （SST-DES）	/	0.31~0.35	0.51
Rockwell（势流）	/	0.38	/
Ahuja & Mendoza^[22]（理论）	/	0.35	0.505
Howe^[23]（理论）		0.39

湍流模型	$ S {t}_{1,\text{Sub}} $	$ S {t}_{1} $	$ S {t}_{2} $	$ S {t}_{2,\text{Mul}} $
DFSEM-WMHRL（NUWA）	0.187	0.375	0.54	1.07
SST-IDDES（STAR-CCM+）	0.189	0.379	/	/

湍流模型	$ S {t}_{1,\text{Sub}} $	$ S {t}_{1} $	$ S {t}_{2} $	$ S {t}_{2,\text{Mul}} $
DFSEM-WMHRL（NUWA）	0.187	0.375	0.54	1.07
SST-IDDES（STAR-CCM+）	0.189	0.379	/	/

基于湍流入口壁面模化RANS/LES混合方法的不可压空腔流数值模拟研究

PDF下载

蓝林欣 ¹ , 韩盼盼 ¹ , 季梦 ¹ , 尤云祥 ¹^,² , 邱小平 ³ , 马乔 ³ , 吴凯健 ³

船舶力学 | 流体力学 2026,30(2): 218-234

收起

船舶力学 |流体力学 2026 , 30 (2) : 218 -234

基于湍流入口壁面模化RANS/LES混合方法的不可压空腔流数值模拟研究

全屏

蓝林欣¹, 韩盼盼¹, 季梦¹, 尤云祥¹^,², 邱小平³, 马乔³, 吴凯健³

作者信息

^1.上海交通大学海洋工程国家重点实验室，上海 200240

^2.上海交通大学三亚崖州湾深海科技研究院，海南三亚 572000

^3.上海君昱信息科技有限公司，上海 200240

通讯作者:

韩盼盼（1990–），男，博士，助理研究员，通讯作者，E-mail：hanpanpan@sjtu.edu.cn

作者简介:

蓝林欣（1993–），男，博士研究生

尤云祥（1963–），男，教授，博士生导师

Numerical investigation on incompressible cavity flows with turbulent incoming flow based on a wall-modeled RANS/LES method

Lin-xin LAN¹, Pan-pan HAN¹, Meng JI¹, Yun-xiang YOU¹^,², Xiao-ping QIU³, Qiao MA³, Kai-jian WU³

Affiliations

^1.State Key Laboratory of Ocean Engineering, Shanghai Jiao Tong University, Shanghai 200240, China

^2.Yazhou Bay Institute of Deepsea Sci-tech, Shanghai Jiao Tong University, Sanya 572000, China

^3.Shanghai Junyu Information Technology Limited Co., Shanghai 200240, China

出版时间: 2026-02-15 doi: 10.3969/j.issn.1007-7294.2026.02.004

文章导航

摘要

收起

空腔流动广泛存在于水下航行器中，如水下航行器的流水孔就属于典型的空腔流动。在空腔流中，存在剪切层K-H不稳定性和腔内回流耦合作用等复杂湍流现象，目前利用CFD来研究此类现象仍面临诸多困难。首先，由于空腔上游来流通常为完全湍流状态，因此湍流入流条件的准确设置是高置信度地解析空腔内湍流相干结构的关键之一。其次，空腔流中除了存在一阶振荡模态的频率成分外，还存在二阶振荡模态的频率成份，而后者迄今仍是CFD数值模拟的难点。有鉴于此，本文发展了一种基于零散度湍流入口模型（简称DFSEM）的空腔流壁面模化RANS/LES混合（简称WMHRL）方法，将其简称为DFSEM-WMHRL模型。通过对$ {Re}_{\tau }=395 $下的槽道流的系列数值模拟研究，获得了能高置信度解析槽道流二阶统计量的湍动能解析度指标$ {r}_{{\mathrm{k}}} $及RANS/LES转换边界位置的组合条件。在此基础上，对$ {Re}_{{\mathrm{D}}}=3360 $及$ L/D=2 $下的空腔流问题采用DFSEM-WMHRL模型进行了系列数值模拟研究。结果表明，该DFSEM-WMHRL模型不仅能够准确解析空腔流动的二阶统计量特征，还能准确获取空腔自持振荡引起的一阶和二阶精细频谱结构。

关键词

空腔流 / 湍流入口 / 相干结构 / 自持振荡频谱 / 壁面模化大涡模拟

Abstract

收起

Cavity flows are common phenomena for underwater vehicles. For example the drain holes of the submarine are one of the typical cavity flows. There are some complex turbulent phenomena in the cavity, such as shear layer K-H instability, coupling between shear layer and cavity recirculation. So far, many difficulties remain for CFD numerical simulation. Firstly, the setup of an accurate inflow condition is crucial for analyzing turbulent coherent structures inside the cavity with the high fidelity due to the fact that the flow upstream of the cavity is usually in an utterly turbulent state. Secondly, there are frequency components with second-order oscillation modes besides the ones with first-order oscillation modes in the cavity flow. However, the CFD numerical simulation still faces challenges for such second-order oscillation flows. In this paper, a DFSEM-WMHRL method is developed to perform high fidelity CFD simulation for such complex cavity flows, where the WMHRL is a hybrid RANS/LES method with the wall-modelled capability, and the DFSEM is a divergence free synthetic eddy method with the turbulent inlet generation capability. Through a series of numerical simulations on the channel flow with $ {Re}_{\tau }=395 $, the turbulence kinetic energy resolution index r_k and the combined conditions for the RANS/LES hybrid boundary positions are firstly confirmed, indicating the capability of the proposed method for analyzing the second-order statistics of the channel flow with high fidelity. Moreover, numerical simulations were carried out for the cavity flow with $ {Re}_{{\mathrm{D}}}=3360 $ and $ L/D=2 $ based on the proposed method. The results show that the DFSEM-WMHRL method can accurately analyze both second-order statistical characteristics and the fine spectral structures caused by self-sustaining oscillations in the cavity flow.

Key words

cavity flows / turbulent inlet / coherent structures / self-sustained oscillation frequencies / wall-modeled LES simulation

引用本文

蓝林欣, 韩盼盼, 季梦, 尤云祥, 邱小平, 马乔, 吴凯健. 基于湍流入口壁面模化RANS/LES混合方法的不可压空腔流数值模拟研究. 船舶力学, 2026 , 30 (2) : 218 -234 . DOI: 10.3969/j.issn.1007-7294.2026.02.004

Lin-xin LAN, Pan-pan HAN, Meng JI, Yun-xiang YOU, Xiao-ping QIU, Qiao MA, Kai-jian WU. Numerical investigation on incompressible cavity flows with turbulent incoming flow based on a wall-modeled RANS/LES method[J]. Journal of Ship Mechanics, 2026 , 30 (2) : 218 -234 . DOI: 10.3969/j.issn.1007-7294.2026.02.004

正文

收起

0　引　言

收起

空腔流动现象在实际工程应用中非常普遍，如飞机起落架、高速列车接缝、汽车天窗及潜艇流水孔等。由于壁面的缺失，当流体流经空腔时，会在空腔前缘发生边界层分离并形成振荡特性的剪切流动，并在一定条件下产生稳定的自持振荡^[1]，从而引起噪声辐射和结构疲劳。所谓自持振荡，是指振荡仅由入流剪切层不稳定性及腔内流动反馈机制共同作用的结果，即在没有其他外部控制或干扰下也会持续产生。近几十年来，国内外学者针对空腔流动问题开展了大量的实验、理论模型及数值模拟研究工作。

早期针对空腔流动的研究大多集中在高马赫数的可压缩气动领域。大量研究表明，在高马赫数下，空腔的振荡特性主要由空腔的声模态特性主导。Rossiter^[2]基于风洞实验研究了亚声速和跨声速状态下的空腔流动，并首次提出了空腔流动的声反馈机制模型，该模型可以较为准确地预测高马赫数下的空腔振荡频率特性。然而，大量研究表明，对于水等不可压缩流体中的空腔流动，声波的波长一般均大于空腔开口长度。以频率

$ f=100\;{\mathrm{Hz}} $

的声波为例，其水下波长

$ \lambda $

约为14.8 m，远大于潜航器上的空腔开口长度。因此，Rossiter^[2]提出的声反馈机制并不适用于水下空腔的振荡频率预报。

对于不可压缩状态下的空腔流动特征，早期主要通过实验手段进行研究，Ethembabaoglu^[3]通过实验研究了

$ 0.25 \lt L/D \lt 6 $

下的不可压缩空腔流动特性，其中L为空腔长度，D为空腔深度。Ethembabaoglu的结果表明，对于不可压缩空腔流，主要存在由自持振荡引起的1~3阶主振荡模态。Lin等^[6]基于线性无粘理论预测了二维不可压缩空腔的振荡频率，与Ethembabaoglu^[3]的实验结果吻合良好。Pereira等^[4−5]通过实验研究了

$ {Re}_{{\mathrm{D}}}=3360、L/D=2 $

下湍流流经浅空腔的流动，分析了直角状、圆弧状及鼻翼状三类空腔后缘形状对腔内流动的影响，结果表明后缘为鼻状型的空腔垂向应力幅值明显减小，但三类空腔下自持振荡频率基本一致；特别地，针对直角状空腔，上述实验还提供了近乎完整的一阶、二阶平均量及频谱特征等参数。

Lin等^[6]采用高精度PIV方法，研究了

$ {Re}_\text{D}\approx 27\;800 $

下的二维空腔问题，并结合流场瞬时速度及压力变化情况，将空腔的自持振荡过程分为以下四个阶段（如图1所示）：第一阶段，空腔前缘来流与腔内流体形成分离剪切层，该剪切层与空腔后缘碰撞后产生压力反馈；第二阶段，这一反馈向上游传播，并对前缘角处形成的剪切层产生影响；第三阶段，到达前缘的扰动转换为分离剪切层中的波动；第四阶段，该波动在分离剪切层向下游传播过程中持续放大。

随着计算流体力学CFD（Computational Fluid Dynamics）技术的发展，数值模拟逐渐成为研究空腔流特性的重要手段之一。由于雷诺平均方法RANS（Reynolds-Averaged Navier-Stokes）只能预测空腔流场的时均量信息^[7]，不能准确捕捉空腔内流体振荡及湍流脉动等高阶瞬态特征量信息，因此并不适用于空腔流这类复杂流动，而直接数值模拟DNS（Direct Numerical Simulation）和大涡模拟方法LES（Large Eddy Simulation）则可以弥补上述缺陷。DNS需要解析边界层内所有尺度的湍流，对网格分辨率的要求非常高，网格量大，不适合高雷诺数钝体绕流计算。LES通过滤波函数将大尺度涡和小尺度涡分离开，对大尺度涡采用直接模拟，小尺度涡通过模型进行封闭，与DNS相比网格量要少很多。Chang等^[8]采用LES方法对

$ {Re}_\text{D}=3360、L/D=2 $

下的空腔流动进行了数值模拟计算，其湍流来流工况在流向平均速度、雷诺应力及频谱特征上均与Pereira等^[4]的实验结果吻合良好。特别地，Chang等在空腔频谱的结果中，还发现了一个低频特征，该频率与一阶自持振荡的次谐波频率相近。

然而，不论是DNS方法或是LES方法，为了准确求解空腔湍流特征，两者对网格数量及数值求解精度的要求都很高。以Chang等^[8]采用的LES数值模拟为例，为了准确地解析空腔内的湍流特征，全域网格均需以Kolmogorov 尺度为基准，因此其对应的

$ {Re}_\text{D}=3360 $

下的空腔流计算网格量达到了1376万。

在实际工程应用上，对水下潜艇等大型水下潜航器，其表面空腔的长度雷诺数

$ {Re}_\text{L} $

一般为10⁴~10⁶。结合Spalart^[10]对10⁵~10⁶下钝体绕流不同湍流模型需要的网格量估计，以

$ {Re}_\text{L}={10}^{6} $

的空腔流动为例， LES所需网格的数量级将达到10¹¹量级。显然，这样的网格量对于实际工程应用计算而言是不适用的。为了解决LES计算量巨大的问题，在过去的二十余年中，各种RANS/LES混合方法相继被提出，在很大程度上解决了LES计算成本高昂的问题，其中IDDES方法在工程领域的应用最为广泛。

IDDES方法通过引入混合函数

$ \tilde{{f}_{d}} $

来实现从RANS到LES的转换，但一般情况下其具有80%湍动能解析能力的LES模式需要在边界层外才能被激活。然而，空腔流场中入流剪切层的K-H不稳定性和腔内回流与剪切层相互作用的复杂湍流特征均发生在边界层内，因此IDDES方法对这类复杂流动问题一般不适用。

本团队成员季梦等^[12]在研究亚临界雷诺数Re=3900下的圆柱绕流问题时，通过修改混合函数

$ \tilde{{f}_{d}} $

，发展了一种壁面模化RANS/LES混合方法WM-HRL（Wall Modeled Hybrid RANS/LES），该方法与IDDES等传统RANS/LES混合方法的不同之处在于，可通过调节网格特征尺度灵活地变换RANS/LES转换区域的位置。对K-H不稳定性发生区域，可通过调整网格使其位于至少具有80%湍动能解析能力的完全LES区，从而实现对钝体绕流的复杂湍流信息进行高置信度解析，但该方法仅适用于来流为层流状态的情况，而对空腔流问题，其来流通常为湍流状态。

Chang等^[8]采用LES方法研究浅空腔问题时，对比了入流为层流和湍流两种情况下的浅空腔流动性，发现在入流为层流时，流体流经空腔前缘产生流动分离后，空腔内的流动仍以层流为主，空腔内的湍流应力水平也仅为实验状态下的20%。此外，在入流为层流来流时，空腔自持振荡引起的特征频率明显失真。由此可见，在入流为层流来流时，即便是高解析度的LES方法也不能模拟实际的空腔流动特性。因此，对空腔流问题，其入流必须采用湍流来流进行模拟。

针对入流为湍流来流的问题，Jarrin等^[13]提出了一种合成涡方法SEM（Synthetic Eddy Method），该方法在入口平面引入人工涡，通过涡的强度来表征速度脉动，并以

$ {Re}_{\tau }=395 $

下的槽道湍流为对象，获得了高置信度的槽道湍流特性。SEM方法也是当前CFD商业软件STAR-CCM+中采用的主要湍流来流生成方法。然而，SEM方法虽然可以生成湍流来流，但仍需要较长的入流段才能形成充分发展的湍流。为克服这个缺陷，Poletto等^[14]提出了一种零散度合成涡方法DFSEM（Divergence Free Synthetic Eddy Method），该方法在SEM方法的基础上，加入了零散度限制，可大大缩短入流段的距离。

有鉴于此，本文将季梦等^[12]发展的WM-HRL模型与Poletto等^[14]提出的DFSEM模型进行结合，发展成一种能够处理入流为湍流的空腔流问题的壁面模化RANS/LES混合方法，简称为DFSEM-WMHRL模型，并以

$ {Re}_\text{D}=3360、L/D=2 $

下完全湍流流经空腔的流动特性问题为研究对象，采用该模型进行系列数值模拟与分析研究，通过系统地分析空腔流场的一阶和二阶统计量、自持振荡频率以及空腔内平均流线等定量与定性特性，对该CFD模型捕捉空腔流动精细特性的可行性和有效性进行评估研究。

1　理论模型

收起

考虑不可压缩流体下的空腔流动问题。为此，建立直角坐标系为

$ O\text{-}{x}_{1}{x}_{2}{x}_{3} $

，其中

$ O{x}_{1} $

轴与来流方向一致，

$ O{x}_{3} $

轴垂直向上为正，

$ u=({u}_{1},{u}_{2},{u}_{3}) $

为流体运动的速度矢量。由于流体是不可压缩的，在RANS框架下，空腔流场的控制方程可写成如下

$ \frac{\partial {\overline{u}}_{k}}{\partial {x}_{k}}=0 $

(1)

$ \frac{\partial {\overline{u}}_{i}}{\partial t}+{\overline{u}}_{k}\frac{\partial {\overline{u}}_{i}}{\partial {x}_{k}}=-\frac{1}{{\rho }_{0}}\frac{\partial \overline{p}}{\partial {x}_{i}}+\frac{\partial }{\partial {x}_{k}}\left(\nu \frac{\partial {\overline{u}}_{i}}{\partial {x}_{k}}\right)\text+\frac{\partial {\tau }_{\text{k}i}}{\partial {x}_{k}} $

(2)

式中：变量上面的横线符号“—”表示雷诺时均，

$ {\rho }_{0} $

为流体密度，

$ p $

为流体的扰动压力，

$ \nu $

为流体的运动粘性系数，

$ {\tau }_{ki} $

为雷诺应力，其在方程（1）~（2）中是未知的，采用如下Boussinesq近似进行计算。

$ {\tau }_{ki}=2{\nu }_{{\mathrm{t}}}{\overline{S}}_{ki}-\frac{2}{3}k{\delta }_{ki} $

(3)

式中：

$ k=\varepsilon /\omega $

为湍动能，

$ \omega $

为比耗率，

$ \varepsilon $

为耗散率，

$ {\overline{S}}_{ki} $

为形变率张量。

本文采用SST k–ω模型^[15]对方程（1）~（3）进行封闭，具体形式如下

$ \frac{\partial k}{\partial t}+{\overline{u}}_{i}\frac{\partial k}{\partial {x}_{i}}={P}_{k}-{\beta }^{*}k\omega +\frac{\partial }{\partial {x}_{i}}\left[\left(\nu \text+{\sigma }_{k}{\nu }_{{\mathrm{t}}}\right)\frac{\partial k}{\partial {x}_{i}}\right] $

(4)

$ \frac{\partial \omega }{\partial t}+{\overline{u}}_{i}\frac{\partial \omega }{\partial {x}_{i}}=\frac{\gamma }{{\nu }_\text{t}}{P}_{k}-\beta {\omega }^{2}+\frac{\partial }{\partial {x}_{i}}\left[\left(\nu \text+{\sigma }_{\omega }{\nu }_\text{t}\right)\frac{\partial \omega }{\partial {x}_{i}}\right]+2\left(1\text-{F}_{1}\right)\frac{{\sigma }_{\omega 2}}{\omega }\frac{\partial k}{\partial {x}_{i}}\frac{\partial \omega }{\partial {x}_{i}} $

(5)

式中：

$ {P}_{k} $

为产生项，

$ {\nu }_\text{t} $

为湍流涡粘系数，其他系数见文献[15]。

涡粘系数

$ {\nu }_\text{t} $

由下式给出

$ {\nu }_\text{t}=\frac{{\alpha }_{1}k}{\max \left({\alpha }_{1}\omega,\left| \overline{S}\right| {F}_{2}\right)} $

(6)

式中：

$ {\alpha }_{1} $

为模型系数，取值为0.31；F₂为混合函数，

$ \left| \overline{S}\right| $

为形变率张量的幅值。

对SST k–ω模型，其IDDES^[11]模型通过修改方程（4）中的色散项而建立

$ \frac{\partial k}{\partial t}+{\overline{u}}_{i}\frac{\partial k}{\partial {x}_{i}}={P}_{k}-\frac{{k}^{3/2}}{{\tilde{l}}_{{\mathrm{hyb}}}}+\frac{\partial }{\partial {x}_{i}}\left[\left(v+{\sigma }_{k}{v}_\text{t}\right)\frac{\partial k}{\partial {x}_{i}}\right] $

(7)

式中：

$ {\tilde{l}}_{{\mathrm{hyb}}} $

为IDDES的混合长度尺度，具体形式如下

$ {\tilde{l}}_{{\mathrm{hyb}}}={\tilde{f}}_{{\mathrm{hyd}}}\cdot \left(1+{f}_{e}\right)\cdot {l}_{{\mathrm{RANS}}}+\left(1-{\tilde{f}}_{{\mathrm{hyd}}}\right)\cdot {l}_{{\mathrm{LES}}} $

(8)

Gritskevich等^[16]指出，上式中引入提升函数

$ {f}_{e} $

旨在增大湍流涡粘系数的值，但这种人为增大湍流涡粘系数的方式并不具有普适性，因此建议取消该提升函数，即将式（8）修改为

$ {l}_{{\mathrm{IDDES}}}={\tilde{f}}_{{\mathrm{hyd}}}\cdot {l}_{{\mathrm{RANS}}}+\left(1-{\tilde{f}}_{{\mathrm{hyd}}}\right)\cdot {l}_{{\mathrm{LES}}} $

(9)

式中：

$ {l}_{{\mathrm{RANS}}} $

为RANS尺度，

$ {l}_{{\mathrm{LES}}} $

为LES尺度，它们可分别表示为

$ {l}_{{\mathrm{RANS}}}={k}^{1/2}/\left({\beta }^{\text{*}}\omega \right),\quad {l}_{{\mathrm{LES}}}={C}_{{\mathrm{DES}}}\Delta $

(10)

式（9）中，

$ {\tilde{f}}_{{\mathrm{hyd}}} $

为混合函数，在当前所见的各类文献（包括商业软件STAR-CCM+等）中采用的均是由Shur等^[11]提出的混合函数

$ {\tilde{f}}_{d} $

的形式，其具体表达式如下

$ {\tilde{f}}_{{\mathrm{hyd}}}=\max \left\{{f}_{B}\text{，}{f}_{dt}\right\} $

(11)

$ {f}_{B}=\min \left(2\exp \left(-9{\alpha }^{2}\right),1.0\right),\alpha =0.25-\frac{{d}_{{\mathrm{w}}}}{{h}_{\max }} $

(12)

$ {f}_{dt}=\tanh \left[{\left(c{d}_{1}\cdot {r}_{dt}\right)}^{c{{d}_{2}}}\right]\text{，}{r}_{dt}=\frac{\nu \text+{v}_\text{t}}{{\kappa }^{2}d_{\text{w}}^{2}\sqrt{\left(\partial {\overline{u}}_{k}/\partial {x}_{i}\right)\left(\partial {\overline{u}}_{k}/\partial {x}_{i}\right)}} $

(13)

式中：

$ {d}_\text{w} $

为计算点到壁面的距离，

$ {h}_{\max } $

为计算单元网格的最大尺度，冯卡门常数

$ \kappa =0.41 $

，

$ c{d}_{1} $

取值为20，

$ c{d}_{2} $

取值为3。

对如上构造的混合函数

$ {\tilde{f}}_{{\mathrm{hyd}}} $

，在通常情况下，其具有80%湍动能解析能力的LES模式需要在边界层外才能被激活。季梦等^[12]指出，对亚临界雷诺数圆柱绕流的K-H不稳定性问题，由于其剪切层K-H不稳定性结果发生在对数律层区内，因此IDDES难以高置信度地解析这类复杂湍流特征，于是提出了一个仅与当地网格空间分布尺寸相关的湍动能解析度指标

$ {r}_{k} $

，并在此基础上提出了一种新的混合函数

$ {\tilde{f}}_{{\mathrm{hyd}}} $

，相应地发展了一种WM-HRL方法。其混合函数

$ {\tilde{f}}_{{\mathrm{hyd}}} $

的具体形式如下

$ {r}_{k}=\min \left[1.0,{\left(0.4\frac{{h}_{\max }}{{d}_\text{w}}\right)}^{2/3}\right] $

(14)

$ {\tilde{f}}_{{\mathrm{hyd}}}=\text{1}-{\mathrm{sign}}\left[{r}_{k\text{1}}-{m}_{k}\right]{f}_{s},\,\,\,\,\,\,\,\,{\mathrm{m}}_{k}=\text{min}\left[{r}_{k\text{1}},{r}_{k}\right] $

(15)

$ {f}_{s}({m}_{k})=\text{min}\left[\text{2exp}\left(-11{\alpha }^{\text{2}}\right),\text{1.0}\right] $

(16)

$ \alpha \text{=min}\left(-\text{0.25,}\,\,\,\,\,\frac{\text{1}}{\text{2}({r}_{k\text{2}}-{r}_{k\text{1}})}{m}_{k}+\frac{{r}_{k\text{1}}-\text{3}{r}_{k\text{2}}}{\text{4}({r}_{k\text{2}}-{r}_{k\text{1}})}\right) $

(17)

式中：

$ {r}_\text{k1} \gt {r}_\text{k2} $

为自定义参数，且取值在0和1之间。

对一般的钝体绕流问题，其入流可设置为层流来流条件。此时，通过合理地设置

$ {r}_\text{k1} $

及

$ {r}_\text{k2} $

的值，一方面可保证RANS区域被充分遮蔽，另一方面又可保证进入完全LES区域时，具有足够的湍动能解析度用以解析复杂湍流特征。

然而，对空腔流动问题，空腔前缘的入流一般为湍流来流条件，为了采用上述WM-HRL方法来研究空腔流问题，还需要解决入流为湍流来流的设置问题。合成湍流方法是实现入流为湍流来流的主要手段之一，其基本思想是将空间速度分解为时间平均速度

$ {\overline{u}}_{i} $

和脉动速度

$ u_{i}^{\prime} $

两部分，如下式所示

$ {u}_{i}={\overline{u}}_{i}+u_{i}^{\prime} $

(18)

Jarrin等^[13]提出了合成涡法SEM（Synthetic Eddy Metod），其中对脉动速度

$ u_{i}^{\prime} $

进行如下分解

$ u_{i}^{\prime}\left(x\right)=\frac{1}{\sqrt{N}}\sum\limits_{k=1}^{N}{a}_{ij}\varepsilon _{j}^{k}f_{\sigma }^{k}\left(\frac{x-{x}^{k}}{{\sigma }^{k}}\right) $

(19)

式中：N为涡的数量，

$ {x}^{k} $

为第k个涡的中心位置，

$ {\sigma }^{k} $

为涡中心的湍流长度尺度，

$ f_{\sigma }^{k} $

为适当的形函数，

$ \varepsilon _{j}^{k} $

是均值为零的随机数表征涡强度，且

$ \left\langle \varepsilon _{j}^{k}\varepsilon _{j}^{k}\right\rangle =1 $

。

$ {a}_{ij} $

为系数矩阵，表达式为

$ {a}_{ij}=\left[\begin{matrix}\sqrt{{R}_{11}} & 0 & 0\\{R}_{21}/{a}_{11} & \sqrt{{R}_{22}-a_{21}^{2}} & 0\\R{}_{31}/{a}_{11} & \left({R}_{32}-{a}_{22}{a}_{31}\right)/{a}_{22} & \sqrt{{R}_{33}-a_{31}^{2}-a_{32}^{2}}\end{matrix}\right] $

(20)

式中：

$ {R}_{ij} $

为雷诺应力张量单元。

Poletto等^[14]指出，采用SEM方法生成入口湍流时，仍需要较长的进流段长度

$ {L}_{{\mathrm{in}}} $

才能获得充分发展的入口湍流场，这意味着仍需要耗费大量网格和计算资源用于进流段计算。有鉴于此，Poletto等在SEM的基础上，提出了一种零散度合成涡方法DFSEM（Divergence Free Synthetic Eddy Method），其实现方法是在散度为零的条件下，将原始SEM方法应用于涡量场，再通过求解涡量场的旋度得到速度场，即

$ \nabla \times {\omega }^{\prime}=\nabla \left(\nabla \cdot {u}^{\prime}\right)-{\nabla }^{2}{u}^{\prime} $

(21)

对不可压缩流体，方程（21）的右端第一项为零，因此式（21）可改写为

$ \nabla \times {\omega }^{\prime}=-{\nabla }^{2}{u}^{\prime} $

(22)

进一步地，脉动速度场可写成^[17]

$ u_{i}^{\prime}\left(x\right)=\frac{1}{\sqrt{N}}\sum\limits_{k=1}^{N}\frac{{q}_{\sigma }\left(\left| {r}^{k}\right| \right)}{{\left| {r}^{k}\right| }^{3}}{r}^{k}\times {\alpha }^{k} $

(23)

式中：

$ {r}^{k}=\left(x-{x}^{k}\right)/{\sigma }^{k} $

，

$ \alpha _{i}^{k} $

是平均值为零的随机数，表示涡的强度。

$ {q}_{\sigma }\left({\left| r\right| }^{k}\right) $

为形函数，表达式如下^[18]

$ {q}_{i}=\begin{cases} {\sigma }_{i}\left[1-{\left({d}^{k}\right)}^{2}\right], & {d}^{k} \lt 1\\0, & {\mathrm{otherwise}}\end{cases} $

(24)

式中：

$ {d}^{k}=\sqrt{{\left(r_{j}^{k}\right)}^{2}} $

，

$ {\sigma }_{i} $

为各向异性的长度尺度。

进一步地，Poletto等表示，对于密度均匀流体，

$ {\sigma }_{i} $

可采用平均长度尺度

$ {\sigma }_{{\mathrm{avg}}} $

代替，表达式如下

$ {\sigma }_{{\mathrm{avg}}}=\min \left({k}^{3/2}/\varepsilon,\kappa \delta,\max \left(\Delta x,\Delta y,\Delta z\right)\right) $

(25)

式中：

$ {k}^{3/2}/\varepsilon $

为局部长度尺度，一般可由RANS方法获得，

$ \delta $

为槽道半高，

$ \kappa =0.40 $

为冯卡门常数。

至此，利用如上所述DFSEM方法，在给定来流平均速度分布

$ {\overline{u}}_{i} $

、雷诺应力分布

$ {R}_{ij} $

以及相应的湍流长度尺度分布

$ {\sigma }_{{\mathrm{avg}}} $

后，即可在入口处生成完全湍流来流。

本文利用作者团队自研CFD软件NUWA：FLOW_UV，在该软件中已有WM-HRL模型的基础上，新植入了DFSEM湍流入口模型，并开发了DFSEM-WMHRL模型的计算模块。该自研CFD软件系统采用有限体积法（FVM）求解控制方程，压力速度求解采用PISO算法，对流项及扩散项的空间离散采用二阶中心差分格式，时间离散格式为二阶隐式格式。

2　结果与分析

收起

为验证本文新开发的DFSEM-WMHRL模型的有效性，首先需要验证利用该模型来模拟入流为湍流来流的有效性。有鉴于此，在本文2.1节中，首先以

$ {Re}_{\tau }=395 $

下的槽道流为对象，进行数值模拟与分析研究。在此基础上，本文2.2节中，以

$ {Re}_\text{D}=3360、L/D=2 $

下完全湍流流经空腔的流动特性问题为研究对象，系统地分析与评估该模型捕捉空腔流动精细特性的可行性和有效性。

在下列数值模拟中，边界条件设置如下：入口采用DFSEM方法生成合成湍流，即

$ u=\overline{u}+{u}^{\prime} $

，湍流入口信息参考Moser等^[20] 基于DNS的数值模拟结果。出口设置为零压力出口，即

$ \overline{p}=0 $

。初始湍动能按湍流强度

$ I=0.2\% $

设定，即

$ k=3{\left({U}_{\mathrm{\infty }}I\right)}^{2}/2 $

。比耗率按

$ \omega =k/{\nu }_\text{t} $

设定，其中初始涡粘比为

$ {\nu }_\text{t}/\nu =0.1 $

。计算区域的上下两个边界设置为无滑移壁面，即

$ {\overline{u}}_{x}={\overline{u}}_{y}={\overline{u}}_{z}=0 $

，湍动能设置为

$ k=0 $

，比耗率

$ \omega $

设置为

$ \omega =6\nu / {\beta }^{*}d_{{\mathrm{w}}}^{2} $

。计算区域的两个垂直侧面设置为对称边界。

2.1　槽道流

对空腔流动问题，其进流段与槽道流动相似。为此，本小节首先以

$ {Re}_{\tau }=395 $

下的槽道流为对象，采用DFSEM-WMHRL模型，对其湍流流动特性进行数值模拟与分析研究。计算域为如图2所示的矩形区域，其中槽道高度为2δ，流向长度为12δ，展向宽度为3δ，流动方向为

$ Ox $

轴正向，y为垂向，z为展向。

边界条件设置如下：在入口边界，速度入口设置为DFSEM合成湍流入口，湍动能按湍流强度

$ {I}=0.2\mathrm{\% } $

设定，即

$ k=3{\left({U}_{\mathrm{\infty }}I\right)}^{2}/2 $

；比耗率按

$ \omega =k/{v}_\text{t} $

设定，其中

$ {v}_\text{t}/v=10.0 $

。在出口边界，设置为零压力出口，即

$ \overline{p}=0 $

。两个垂直侧面设置为对称边界条件。进流段与出流段上下边界均设置为无滑移条件，即

$ {\overline{u}}_{x}={\overline{u}}_{y}={\overline{u}}_{z}=0 $

，湍动能设置为

$ k=0 $

，比耗率

$ \omega $

设置为

$ \omega =6v/{\beta }^{*}d_{{\mathrm{w}}}^{2} $

。

对DFSEM-WMHRL模型的RANS模式，其边界层内第一层网格的无量纲壁面距离

$ {y}^{+} $

应满足

$ {y}^{+}\leq 1 $

。其中，

$ {y}^{+} $

定义为

$ {y}^{+}=y\sqrt{{\tau }_{{\mathrm{w}}}}/\left(\nu \sqrt{\rho }\right) $

，

$ y $

为网格单元中心与壁面的距离，

$ {\tau }_{{\mathrm{w}}} $

为壁面切应力。在本文中，第一层网格的无量纲壁面距离取为0.5，网格法向增长率取为1.16。在此基础上，采用ANSYS ICEM软件对槽道进行网格划分，单元网格为六面体，其纵剖面网格划分如图2（b）所示。在数值计算中，无量纲化时间步长

$ \Delta {t}^{*}(=\Delta t{U}_{\mathrm{\infty }}/D) $

取值为

$ 2.2\times {10}^{-2} $

，最大库朗数

$ C{o}_{\max }\approx 0.5 $

，统计平均开启时刻为

$ T={70D/U}_{\mathrm{\infty }} $

，统计平均时长约为

$ 200D/{U}_{\mathrm{\infty }} $

，并将结果与DNS方法^[20]及PANS（Partially Averaged Navier–Stokes）方法^[19]的数值结果进行比较分析。

对

$ {Re}_{\tau }=395 $

下的槽道流，其粘性底层的

$ {y}^{+}\leq 10.0 $

，过渡层的

$ {y}^{+} $

为10~40，对数率层的

$ {y}^{+} $

为40~400^[20]。对该槽道流问题，DNS方法^[20]的数值结果表明，其湍流脉动峰值出现在

$ y_{\delta \max }^{+}=14.5 $

附近，位于过渡层内。有鉴于此，在本文中，对DFSEM-WMHRL模型所需要确定的LES边界位置，设置为三种情况，分别为

$ y_{{\mathrm{LES}}}^{+} $

=12、14和16，这三种情况的LES边界位置均位于过渡层内，其中第一和第三种情况的LES边界离

$ {y}^{+}=14.5 $

的边界均较远，而第二种情况的LES边界与

$ {y}^{+}=14.5 $

的边界接近。此外，Klapwijk等^[19]采用LES、PANS及IDDES方法对

$ {Re}_{\tau }=395 $

下的槽道完全湍流问题进行了系列数值模拟计算。其中，PANS方法的结果表明，只有当湍动能解析度

$ {r}_{k}\geq 80\mathrm{\% } $

时，才能准确获取槽道的湍流脉动特征。有鉴于此，在本文中将DFSEM-WMHRL模型中的

$ {r}_\text{k2} $

设置为0.2。

在本文中，槽道流计算区域网格主要依据

$ {r}_\text{k2} $

和

$ y_\text{LES}^{+} $

的组合而设计。其中，当

$ {r}_\text{k2}=0.2 $

，而

$ y_\text{LES}^{+} $

=12、14和16时，在这三种组合下的网格量分别为310万、230万和170万。在网格设置完成后，可根据式（14）计算

$ {r}_{k} $

。在此基础上，可设置不同的RANS边界位置及其相应的

$ {r}_\text{k1} $

组合。在本文中，考虑如下三种情况，分别为

$ y_\text{RANS}^{+}=3.2,{r}_\text{k1}=0.56 $

；

$ y_\text{RANS}^{+}=4.6,{r}_\text{k1}=0.44 $

和

$ y_\text{RANS}^{+}=6.3,{r}_\text{k1}=0.35 $

。其中，这三种情况的RANS边界均位于粘性底层内。为后文陈述简便计，将这三种情况分别记为Case A、Case B和Case C。

在图3中，分别给出了当

$ {f}_\text{k2}=0.2 $

，

$ y_\text{LES}^{+} $

分别等于12、14及16时，在Case A、Case B和Case C三种情况下，槽道流的平均速度

$ {\overline{u}}^{+} $

、流向正应力

$ \overline{u}{\overline{u}}^{\prime} $

、垂向正应力

$ \overline{v}{\overline{v}}^{\prime} $

以及水平剪切应力

$ \overline{u}{\overline{v}}^{\prime} $

分布的计算结果。其中，上述为x=10δ处，展向（z方向）进行平均的结果。同时，图中还与DNS^[20]、LES^[19]和IDDES^[19]三种方法的结果进行了比较，其中DNS方法用了1700万网格，LES方法用了950万网格，IDDES方法用了260万网格。

由图可知，对一阶统计量，DNS^[20]、LES^[19]和IDDES^[19]三种方法的结果均吻合一致，而且当完全LES模式启动边界的解析度满足

$ {r}_\text{k2}=0.2 $

时，在三种不同

$ y_\text{LES}^{+} $

位置情况下均可准确地模拟槽道流的一阶速度分布特征。

对于流向正应力

$ \overline{u}{\overline{u}}^{\prime} $

，当

$ {y}^{+} \gt 180 $

时，LES^[19]法的结果与DNS法^[20]的结果出现明显偏差，而当

$ {y}^{+} \gt 260 $

时，IDDES方法^[20]的结果与DNS法^[20]的结果出现明显偏差。对本文的DFSEM-WMHRL模型，除了峰值外，Case A~C三种情况的结果均与DNS法^[20]结果一致，且对Case C（即

$ y_\text{LES}^{+} \lt 14 $

且

$ {r}_\text{k1} \lt 0.5 $

），本文所得峰值结果与DNS法^[20]的结果也吻合良好。其中，

$ y_\text{LES}^{+} \lt 14 $

表示LES法完全启动边界位置所需的无因次壁面距离，需要小于流向正应力峰值所对应的无因次壁面距离（即

$ y_{\delta \max }^{+}\approx 14.5 $

^[20]）。

对水平剪切应力

$ \overline{u}{\overline{v}}^{\prime} $

，LES^[19]和IDDES^[19]法的结果与DNS法^[20]的结果也均出现明显偏差。对本文的DFSEM-WMHRL模型，对Case A~C三种工况所得曲线的变化特征均与DNS^[20]的结果一致，但峰值所对应的

$ {y}^{+} $

值与DNS法^[20]结果相应的

$ {y}^{+} $

值相比，均偏大，其中Case C所得曲线峰值所对应的

$ {y}^{+} $

值与DNS^[20]结果相应的

$ {y}^{+} $

值最为接近，两者的误差为15%。

对于垂向正应力

$ \overline{v}{\overline{v}}^{\prime} $

，LES^[19]和IDDES^[19]两种方法的结果与DNS^[20]的结果均出现明显偏差。对本文的DFSEM-WMHRL模型，Case A和Case B的结果与DNS法^[19]的结果也均出现明显偏差，对于Case C，当

$ {y}^{+} \gt 140 $

时，其结果与DNS法^[20]的结果吻合良好，当

$ {y}^{+} \lt 140 $

时，其结果与DNS法^[20]的结果略有差异，最大误差为9.5%。

下面以Case C为对象，研究完全LES模式开启边界湍动能的解析度

$ {r}_\text{k2} $

对槽道流一阶及二阶统计量的影响，并与DNS法^[20]的结果进行比较，如图4所示。其中，

$ {r}_\text{k2} $

分别取为0.2、 0.25和0.3,即完全LES模式开启边界湍动能的解析度分别为80%、75%和70%。由图可知，三个

$ {r}_\text{k2} $

的取值对一阶统计量及二阶流向正应力计算结果的影响并不明显，但对二阶垂向正应力及水平剪切应力计算结果的影响明显。结果表明，为准确模拟槽道流一阶和二阶统计量特性，DFSEM-WMHRL模型的完全LES模式开启边界湍动能的解析度至少需要达到80%，即

$ {r}_\text{k2}\leq 0.2 $

。80%的湍动能解析度已具备解析复杂湍流相干结构的能力。若想进一步提高边界解析度，则需要细化RANS区及RANS/LES混合区内的网格，使边界层内网格数量激增。对RANS区域，其对湍动能的解析度为零，网格数量的增大并不能提升其对湍流的解析能力，反而造成大量网格资源和计算资源的浪费。

综上可见，对DFSEM-WMHRL模型，为准确获取

$ {Re}_{\tau }=395 $

下的槽道湍流的一阶及二阶统计量特性，需满足以下三个要求：（1）完全LES开启边界位置应位于二阶正应力峰值对应的无量纲壁面距离之前，即

$ y_\text{LES}^{+} \lt y_{{\sigma }_{\max }}^{+} $

；（2）完全LES开启边界的湍动能解析度应不小于80%，即

$ {r}_\text{k2}\leq 0.2 $

；（3） RANS结束边界的湍动能解析度至少需要50%，即

$ {r}_\text{k1}\leq 0.5 $

。

2.2　空腔流

本小节对

$ {Re}_\text{D}=3360 $

下具有完全湍流的入口来流的空腔流问题进行数值模拟研究与分析。其中，

$ Re={U}_{\mathrm{\infty }}D/v $

，

$ {U}_{\mathrm{\infty }} $

为空腔来流平均速度，

$ {D} $

为空腔高度。参考Pereira等^[4]的实验布置及Chang等^[8]的LES数值模拟设置，本文计算域如图5所示，其中，L/D=2，H/D≈6.09，W/D=6，L为空腔长度，H为空腔进流段高度，W为空腔宽度。另外，空腔进流段长度为10D，出流段长度为4D。流动方向为

$ Ox $

轴正向，y为垂向，z为展向。

边界条件设置如下：在入口边界，速度入口设置为DFSEM合成湍流入口，湍动能按湍流强度

$ {I}=0.2\mathrm{\% } $

设定，即

$ k=3{\left({U}_{\mathrm{\infty }}I\right)}^{2}/2 $

；比耗率按

$ \omega =k/{v}_\text{t} $

设定，其中

$ {v}_\text{t}/v=10.0 $

。在出口边界，设置为零压力出口，即

$ \overline{p}=0 $

。两个垂直侧面设置为对称边界条件。进流段与出流段上下边界及空腔底部边界，均设置为无滑移条件，即

$ {\overline{u}}_{x}={\overline{u}}_{y}={\overline{u}}_{z}=0 $

，湍动能设置为

$ k=0 $

，比耗率

$ \omega $

设置为

$ \omega =6v/{\beta }^{*}d_{{\mathrm{w}}}^{2} $

。

对DFSEM-WMHRL模型的RANS模式，其第一层网格的无量纲壁面距离取为0.5，网格法向增长率取为1.1。在此基础上，采用ANSYS ICEM软件对空腔进行网格划分，单元网格为六面体，其纵剖面网格划分如图5所示。在数值计算中，无量纲化时间步长

$ \Delta {t}^{*}(=\Delta t{U}_{\mathrm{\infty }}/D) $

取值为

$ 2.1\times {10}^{-3} $

，最大库朗数

$ C{o}_{\max }\approx 0.3 $

，计算总时长为300D/U，并对后200D/U的流动参数进行统计平均，以获取空腔流场的积分量及流场统计量等信息。

由Chang等^[8]LES数值模拟结果可知，空腔进流段的摩擦雷诺数为

$ {Re}_{\tau }=505 $

。DNS^[21]的结果表明，在该摩擦雷诺数下槽道流二阶正应力峰值所对应的无因次壁面距离

$ y_{{\sigma }_{\max }}^{+} $

约为14.5。根据2.1节的结果，对空腔流问题，本文将DFSEM-WMHRL模型所需的完全LES模式开启边界位置

$ y_\text{LES}^{+} $

取为12.0。此外，根据2.1节的结果，本文将DFSEM-WMHRL模型所需的完全LES模式开启边界的湍动能解析度

$ {r}_\text{k2} $

取为0.2，即保证DFSEM-WMHRL模型在其开启LES模式后对湍动能的解析能力至少达到80%。

在本文中，对空腔流计算区域（如图5所示）的网格，依据如上给定的

$ y_\text{LES}^{+} $

和

$ {r}_\text{k2} $

值进行设置。对空腔流的数值模拟，其进流段和去流段网格依据上述设置可准确地模拟完全湍流来流条件。但对于空腔内部及空腔上方的网格设置仍需进行相关的网格收敛性验证。有鉴于此，本文设计了三套网格方案进行网格收敛性验证。表1为三套网格方案参数及相应的数值计算设置。其中，三套网格方案的进流段及去流段网格划分方式一致。

基于表1中的三套网格方案，利用DFSEM-WMRHL模型对空腔流分别进行数值模拟，并选取两个典型站位

$ x/D=1.0 $

及

$ x/D=1.4 $

下的流向平均速度剖面与Pereira等^[5]的试验及Chang等^[8]LES结果进行比对，以检验网格的收敛性。同时，每套网格方案下，对应的时间步长均满足最大库朗数

$ C{o}_{\max }\approx 0.3 $

。图6为三套网格方案下，两个站位处的无因次流向平均速度计算结果。需要说明的是，所有结果均已进行展向（z方向）平均处理。

由图6可知，在Grid 1粗糙网格方案下，空腔不同站位处流向平均速度与实验^[5]及LES^[8]结果均出现较大偏差。而在Grid 2精细网格及Grid 3超精细网格方案下，均可以准确地获取流向平均速度剖面特性。不难看出，选用Grid 2网格划分方案可以兼顾数值计算精度和效率，因此本文选择Grid 2网格划分方案用于后续的数值模拟计算。

在Grid 2网格划分方案下，本文除了利用DFSEM-WMRHL模型对空腔流进行数值模拟研究外，还利用商业软件STAR-CCM+中的SST-IDDES模型进行了数值模拟研究，以评估该RANS/LES混合模型对空腔流的解析能力。Chang等^[8]利用LES对该空腔流进行了数值模拟研究，所用网格数为1376万。此外，Chang等^[9]还利用SST-DES模型对空腔流进行了数值模拟研究，所用网格数为318万。由此可见，本文所用网格的量约为Chang等^[8]所用的LES网格量的1/4，而与Chang等^[9]所用的SST-DES的网格量相当。

Pereira等^[5]对该空腔流7个不同站位处的流向平均速度（一阶统计量）分布特性进行了试验测量与分析。其中，在空腔入流段设置了1个测点，其位置为

$ x/D=-0.1 $

；在空腔内部设置了5个测点，其位置分别为

$ x/D= $

0.3、0.6、1.0、1.4和1.7；在空腔出流段设置1个测点，其位置为

$ x/D=2.1 $

。Chang等^[8]利用LES对这7个测点处的流向平均速度分布特性进行了数值模拟与分析。

本文分别利用商业软件STAR-CCM+中的SST-IDDES模型，以及NUWA:FLOW_UV中新开发的DFSEM-WMHRL模型对这7个测点处的流向平均速度分布特性进行了数值模拟与分析，同时与Pereira等^[5]的试验结果和Chang等^[8]的LES结果进行了比较分析，结果如图7所示。由图可知，无论是Chang等^[8]的LES模型，还是商业软件STAR-CCM+的SST-IDDES模型，以及NUWA的DFSEM-WMHRL模型，其对空腔7个不同站位处流向平均速度分布特性的数值模拟结果，均与Pereira等^[5]的试验结果吻合良好。

Pereira等^[5]对该空腔流6个不同站位处的二阶统计量分布特性进行了试验测量与分析。其中，在空腔入流段设置了1个测点，其位置为

$ x/D=-0.1 $

；在空腔内部设置了5个测点，其位置分别为

$ x/D= $

0.3、0.6、1.0、1.4和1.7。Chang等^[8]利用LES对这6个测点处的二阶统计量分布特性进行了数值模拟与分析，而Chang等^[9]则利用SST-DES对这6个测点处的二阶统计量分布特性进行了数值模拟与分析。本文分别利用商业软件STAR-CCM+中的SST-IDDES模型，以及NUWA:FLOW_UV中新开发的DFSEM-WMHRL模型对这6个测点处的二阶统计量分布特性进行了数值模拟与分析，同时与Pereira等^[5]的试验结果、Chang等^[8]的LES结果及Chang等^[9]的SST-DES结果进行了比较分析，结果如图8所示。

由图8可知，SST-DES模型的数值结果与试验结果有明显的差异，特别其峰值与试验结果的峰值相差甚大。对SST-IDDES模型，在

$ x/D= $

0.3站位处的峰值与试验结果有较大差异，而在其他站位处的峰值均与试验结果接近，但当

$ y/D \gt 0.2 $

时其数值结果均与试验结果有显著差异。对Chang等^[8]的LES模型和NUWA的DFSEM-WMHRL模型，在6个不同站位处二阶统计量分布的峰值均与试验结果一致，而且当

$ y/D \lt 0.5 $

时两个模型的数值结果均与试验结果吻合，但当

$ y/D\geq 0.5 $

时两个模型的数值结果均与试验结果有一定的差异。

针对本文所述空腔流问题，其流动特征主要包括：（1）流向剪切层K-H不稳定性引起的一阶振荡模态特征及关联的次谐波振荡模态特征；（2）空腔内回流与流向剪切层作用引起的二阶振荡模态特征及关联的高阶振荡模态特征。为研究这两类现象，参考Pereira等^[5]的文献，在空腔内设置了如图10所示的6个测点，用以监测腔内垂向速度的脉动特性，并通过其Lomb谱来获取腔内流体振荡激发的频谱特性。这6个测点的具体坐标信息如表2所示。

针对本文所述空腔流问题，综合作者所能找到的相关文献，总计发现有四类频率成分，如表3所示。其中，

$ S {t}_{1} $

为流向剪切层K-H不稳定性引起的一阶振荡模态无因次频率，

$ S {t}_{2} $

为空腔内回流与流向剪切层作用引起的二阶振荡模态无因次频率，

$ S {t}_{1,{\mathrm{Sub}}} $

为一阶振荡模态频率关联的次谐波无因次频率，而

$ S {t}_{2,{\mathrm{Mul}}} $

为二阶振荡模态频率关联的二倍及三倍无因次频率。其中，所谓无因次频率

$ St $

，定义为

$ St=fD/{U}_{\mathrm{\infty }} $

，f为相应的有因次频率。

由表3可知，Pereira等^[4]的试验只测得1个频率成分

$ S {t}_{2} $

。进一步地，Pereira等^[5]在后续的试验中也只测得2个频率成分

$ S {t}_{2} $

和

$ S {t}_{2,{\mathrm{Mul}}} $

。Chang 等^[9] 利用SST-DES的数值模拟只得到2个频率

$ S {t}_{1} $

和

$ S {t}_{2} $

，而Chang 等^[8] 利用LES的数值模拟也只得到3个频率

$ S {t}_{1,{\mathrm{Sub}}} $

、

$ S {t}_{1} $

和

$ S {t}_{2} $

。此外，利用相关的简化理论模型只能得到1个频率

$ S {t}_{1} $

^[23]。

在本文中，分别利用商业软件STAR-CCM+中的SST-IDDES模型和自研软件NUWA中的DFSEM–WMHRL模型，对该空腔内6个测点处的垂向脉动速度进行了计算，并绘制了相应的Lomb谱，结果如图9所示。由SST-IDDES计算得到的结果可知，在1~3号和5号这4个测点处，只有一个频率成分，即

$ S {t}_{1,{\mathrm{Sub}}} $

。在4号和6号这2个测点处，出现两个频率成分，分别为

$ S {t}_{1,{\mathrm{Sub}}} $

和

$ S {t}_{1} $

。

由DFSEM-WMHRL计算得到的结果可知，在6个测点处，均出现一阶振荡模态频率关联的次谐波无因次频率成分，即

$ S {t}_{1,\text{Sub}} $

。在2号测点处，出现两个频率成分，其中一个为流向剪切层K-H不稳定性引起的一阶振荡模态频率（即

$ S {t}_{1} $

），而另一个为一阶振荡模态频率关联的次谐波频率（即

$ S {t}_{1,\text{Sub}} $

）。在5号测点，则出现了3个频率成分，除了

$ S {t}_{1} $

和

$ S {t}_{1,\text{Sub}} $

外，还出现了空腔内回流与流向剪切层作用引起的二阶振荡模态频率成分（即

$ S {t}_{2} $

）。在3号测点，除了出现

$ S {t}_{1} $

、

$ S {t}_{1,\text{Sub}} $

和

$ S {t}_{2} $

这3个频率成分外，还出现了第四个频率成分，即与二阶振荡模态频率关联的二倍频率成分（即

$ S {t}_{2,\text{Mul}} $

）。

表4中给出了基于SST-IDDES和DFEM-WMHRL模型所得腔内各频率成分无因次值的结果。进一步地，对SST-IDDES模型，其所得两个无因次频率

$ S {t}_{1,{\mathrm{Sub}}} $

和

$ S {t}_{1} $

值，与Chang 等^[8] 利用LES所得相应计算值基本一致，相对误差分别为2.16%和−0.26%。对DFSEM-WMHRL模型，计算共获得四个频率成分，其中前两个无因次频率

$ S {t}_{1} $

和

$ S {t}_{1,\text{Sub}} $

值，也与Chang 等^[8] 利用LES所得相应计算值一致，相对误差分别为1.08%和−1.32%；对二阶振荡频率

$ S {t}_{2} $

，则与Pereira等^[4]的实验结果接近，相对误差为−1.82%；对二阶振荡模态频率关联的二倍频率成分

$ S {t}_{2,\text{Mul}} $

，与Pereira等^[4]的实验结果相比明显偏小，相对误差为−13.71%。

最后，给出空腔流平均流线特征的SST-IDDES和DFSEM-WMHRL模型的数值模拟结果，并与相关文献中的试验^[4]、LES^[8]及SST-DES^[9] 的结果进行比较，结果如图10所示。由图可知，Pereira等^[4]在试验中观测到的两个再循环涡流特征，分别由一个较大的主涡和一个较小的次涡构成，其中主涡为空腔前缘剪切层作用下形成的顺时针方向旋转的大涡，次涡为腔内回流作用下形成的逆时针方向旋转的小涡。Chang等利用LES^[8]及SST-DES^[9]得到的数值结果中均捕捉到主涡及次涡这两个再循环涡流特征，且两个涡流特征的位置与大小均与试验接近。此外，Chang等利用LES^[8]及SST-DES^[9]数值模拟还在空腔后缘底部捕捉到一个极小的再循环涡流特征。本文基于SST-IDDES及DFSEM-WMHRL模型所得的空腔内平均流线特征中，也均准确地捕捉到主涡及次涡特征，且其位置分布和涡结构尺度也与试验^[4]、LES^[8]及SST-DES^[9]结果一致。特别地，SST-IDDES及DFSEM-WMHRL模型也均能在空腔后缘底部捕捉到一个极小的再循环涡流特征。对于试验中未能捕捉到空腔后缘底部存在的极小再循环涡流特征，从图9可以看出，该再循环涡流特征水平位置位于

$ 1.8\leq x/D\leq 2.0 $

，而试验中距离空腔后缘最近的流速测量剖面位于

$ x/D=1.7 $

，因此Pereira等^[4]在试验中无法捕捉到这一特征。

3　结　论

收起

对于空腔流的数值模拟，除了需要准确模拟其相干结构等壁面瞬态湍流场外，还需要准确模拟其进流端的壁面瞬态完全湍流场。前者可以采用季梦等^[12]发展的WM-HRL模型进行数值模拟，但季梦等^[12]所采用的层流入流方式对后者不再适用，需要采用Poletto等^[14]提出的零散度合成涡方法（简称DFSEM）。本文将两者相结合，发展了一种针对具有完全湍流入口条件的壁面模化RANS/LES混合方法（简称DFSEM-WMHRL），开发了相应的CFD计算程序，并已植入作者团队自主开发的NUWA:FLOW_UV软件。

通过对

$ {Re}_{\tau }=395 $

下的槽道完全湍流问题的数值模拟研究与分析，获得了本文所提出的DFSEM-WMHRL模型中完全LES模式启动边界位置

$ y_\text{LES}^{+} $

及其湍动能解析度

$ {r}_\text{k2} $

,以及RANS结束边界位置

$ y_\text{RANS}^{+} $

及其湍动能解析度

$ {r}_\text{k2} $

这四个模型参数取值的准则。在此基础上，对文献[5]中

$ {Re}_\text{D}=3360、L/D=2 $

的不可压空腔流特性进行了数值模拟与分析，得到如下主要结论：

（1）目前常用的SST-DES和SST-IDDES等脱体涡型RANS/LES混合模型在采用相关的合成涡模拟方法后，虽然能够较为准确地模拟空腔一阶统计量的分布特性，但均不能准确模拟空腔的二阶统计量分布特性。对本文所发展的DFSEM-WMHRL模型，在采用相关的合成涡模拟方法后，除了能够准确地模拟空腔一阶统计量的分布特性，还能准确模拟空腔的二阶统计量分布特性。

（2）本文所述空腔流问题存在4个振荡频率成分，目前相关文献中无论是试验研究还是基于LES、SST-DES以及基于势流和简化理论模型的研究，均不能准确获得这4个频率成分的完整且准确的信息，而利用本文所发展的DFSEM-WMHRL模型，可以获得这4个频率成分的完整且准确的信息。

（3）对于空腔流平均流线特征，文献中相关的LES及SST-DES模型以及本文所采用的SST-IDDES及SEM-WMHRL模型，均可以准确地模拟出空腔内三个再循环涡特征，其中再循环涡大小和位置与相关文献中的试验结果一致。

综上可见，本文提出的DFSEM-WMHRL模型可以准确地解析空腔流的湍流相干结构特征及腔内各类自持振荡频率成分特性，为后续高雷诺数下各类浅腔及深腔空腔流的高置信度CFD数值模拟提供了一种新的手段。

参考文献

收起

文献

收起

参考文献引证文献

排序方式：

Rockwell D, Nandascher E. Review-self-sustaining oscillations of flow past cavities[J]. Transactions of the ASME, 1978, 100: 152−165.

Rossiter J E. Wind-tunnel experiments on the flow over rectangular cavities at subsonic and transonic speeds[M/OL]. Aeronautical Research Council. 1964. https://reports.aerade.cranfield.ac.uk.handle.1826.2/4020.

Ethembabaoglu S. On the fluctuationg flow characteristics in the vicinity of gate slots[M]. Division of Hydraulic Engineering, University of Trondheim, Norwegian Institute of Technology, 1973.

Pereira J C F, Sousa J M M. Influence of impingement edge geometry on cavity flow oscillations [J]. AIAA Journal, 1994, 32 (8): 1737−1740.

Pereira J C F, Sousa J M M. Experimental and numerical investigation of flow oscillations in a rectangular cavity[J]. Journal of Fluids Engineering-Transactions of the ASME, 1995, 117(1): 68−74.

Lin J C, Rockwell D. Organized oscillations of initially turbulent flow past a cavity[J]. AIAA Journal, 2001, 39 (6): 1139−1151.

陈浩, 袁先旭, 毕林, 等. 基于RANS/LES混合方法的分离流动模拟[J]. 航空学报, 2020, 41 (8): 123642.

Chen H, Yuan X X, Bi L, et al. Simulation of separated flow based on RANS/LES hybrid method [J]. Acta Aeronautica ET Astronautica Sinica, 2020, 41 (8): 123642. (in Chinese)

Chang K, Constantinescu G, Park S. Analysis of the flow and mass transfer processes for the incompressible flow past an open cavity with a laminar and a fully turbulent incoming boundary layer[J]. Journal of Fluid Mechanics, 2006, 561: 113‒145.

Chang K, Constantinescu G, Park S. Assessment of predictive capabilities of detached eddy simulation to simulate flow and mass transport past open cavities[J]. Transactions of the ASME, 2007, 129: 1372‒1383.

Spalart P R. Strategies for turbulence modelling and simulations[J]. International Journal of Heat & Fluid Flow, 2000, 21(3): 252‒263.

Shur M L, Spalart P R, Strelets M K, et al. A hybrid RANS-LES approach with delayed-DES and wall-modelled LES capabilities[J]. International Journal of Heat & Fluid Flow, 2008, 29(6): 1638‒1649.

季梦, 尤云祥, 韩盼盼, 等. 亚临界区圆柱绕流相干结构壁面模化混合RANS/LES模型[J]. 物理学报, 2024, 73(5): 196‒222.

Ji M, You Y X, Han P P, et al. A wall-modeled hybrid RANS/LES model for flow around a circular cylinder with coherent structures in subcritical Reynolds number regions[J]. Acta Physica Sinica, 2024, 73(5): 196‒222. (in Chinese)

Jarrin N, Prosser R, Uribe J C, et al. Reconstruction of turbulent fluctuations for hybrid RANS/LES simulations using a Synthetic-Eddy Method [J]. International Journal of Heat and Fluid Flow, 2009, 30(3): 435−442.

Poletto R, Craft T, Revell A. A new divergence free synthetic eddy method for the reproduction of inlet flow conditions for LES[J]. Flow, Turbulence and Combustion, 2013, 91(3): 519‒539.

Menter F R. Two-equation eddy-viscosity turbulence modeling for engineering applications[J]. AIAA Journal, 1994, 32: 1598‒1605.

Gritskevich M S, Garbaruk A V, Schütze J, et al. Development of DDES and IDDES formulations for the k-ω shear stress transport model[J]. Flow, Turbulence and Combustion, 2012, 88: 431−449.

Davidson L. Using isotropic synthetic fluctuations as inlet boundary conditions for unsteady simulations[J]. Advances and Application in Fluid Mechanics, 2007, 1(1): 1‒35.

Enciso A, et al. The Biot-Savart operator of a bounded domain[J]. Journal de Mathématiques Pures et Appliquées.

Klapwijk M, et al. Evaluation of scale-resolving simulations for a turbulent channel flow[J]. Computers & Fluids, 2020, 209: 104636.

Moser R D, Kim J, Mansour N N. Direct numerical simulation of turbulent channel flow up to Re_T=590[J]. Physics of Fluids, 1999, 11: 943−945.

Bernardini M, Pirozzoli S, Orlandi P. Velocity statistics in turbulent channel flow up to Re_τ = 4000[J]. Journal of Fluid Mechanics, 2014, 742: 171‒191.

Ahuja K K, Mendoza J. Effects of cavity dimensions, boundary layer, and temperature on cavity noise with emphasis on benchmark data to validate computational aeroacoustic codes: NASA-CR-4653[R]. NASA, 1995.

Howe M S. Edge, cavity and aperture tones at very low mach numbers[J]. Journal of Fluid Mechanics, 1996, 330: 61−84.

2026年第30卷第2期

PDF下载

引用本文

BibTeX

文章信息

doi: 10.3969/j.issn.1007-7294.2026.02.004

接收时间：2024-12-17
首发时间：2026-07-07
出版时间：2026-02-15

补充材料

相关文章

文章信息

作者

出版历史

收稿日期：2024-12-17

基金

作者信息

^1.上海交通大学海洋工程国家重点实验室，上海 200240

^2.上海交通大学三亚崖州湾深海科技研究院，海南三亚 572000

^3.上海君昱信息科技有限公司，上海 200240

通讯作者:

韩盼盼（1990–），男，博士，助理研究员，通讯作者，E-mail：hanpanpan@sjtu.edu.cn

参考文献

分享链接

https://castjournals.cast.org.cn/joweb/cblx/CN/10.3969/j.issn.1007-7294.2026.02.004

分享至

全文二维码

扫描看全文

引用本文

BibTeX

本文的引用情况

2种不同金属材料的力学参数

科 Family	属数 Number of genus	种数 Number of species	占总种数比例 Percentage of total species (%)	属 Genus	种数 Number of species	占总种数比例 Percentage of total species (%)
鹅膏菌科Amanitaceae	2	11	5.26	鹅膏菌属 Amanita	10	4.78
小菇科 Mycenaceae	2	12	5.74	丝盖伞属 Inocybe	5	2.39
多孔菌科 Polyporaceae	8	14	6.70	蜡蘑属 Laccaria	5	2.39
红菇科 Russulaceae	3	23	11.00	小皮伞属 Marasmius	6	2.87
				小菇属 Mycena	11	5.26
				光柄菇属 Pluteus	5	2.39
				红菇属 Russula	17	8.13
				栓菌属 Trametes	5	2.39

关闭全屏

BibTeX
EndNote
RefWorks
TxT