本文摘要：（由ai生成）
本研究通过直接数值模拟（DNS）和耦合的level-set与VOF方法研究了通气空化流。模拟揭示了空泡内的三种特征流动结构，并发现湍动能集中在喷射层。通过湍动能收支分析和谱本征正交分解（SPOD），揭示了不同区域中湍流的产生和耗散机制，以及泄气和涡脱落的相干流结构。本研究深入理解了通气空化流的动力学特性，为相关领域的进一步研究提供了重要参考。
本文通过DNS方法研究通气空化流，使用耦合level-set和VOF方法捕捉气相与水相之间的界面。通气圆盘空化器用于创建空泡，并通过锐界面浸没边界法进行建模。模拟数据全面描述了空化流中的两相流以及泄气和旋涡脱落过程。空气相的平均速度表明存在三种特征流动结构，即剪切层、再循环区和喷射层。湍动能集中在闭合区域的JL，流向湍流波动在SL和JL的横向波动中占主导地位。TKE的收支分析表明，由于速度梯度较高，生成项导致TKE在SL中增加，而在JL中由于流向拉伸效应而减少。在闭合区域中周期性地发生泄气和涡脱落，这两个过程之间的一一对应通过流体速度和体积的谱结果以及空气体积分数的自相关函数得到证实。此外，用谱本征正交分解法分析了相干流结构。我们确定了几种精细的相干结构，包括由开尔文–亥姆霍兹不稳定性引起的SLKH、与大规模旋涡脱落相关的SLout、与小规模旋涡脱落相关的SLin以及与上游湍流对流相关的SLr。本研究补充了前人的研究，详细描述了与空气和水的剧烈混合相关的湍流运动，以及空泡流中不同特征结构之间的复杂相互作用。
近期，美国明尼苏达大学机械工程系和St. Anthony Falls实验室与北京大学工学院湍流与复杂系统国家重点实验室的Han Liu, Zuoli Xiao and Lian Shen等学者利用DNS方法与耦合level-set和VOF方法在2亿网格上利用浸没边界法与两相不可压缩求解器对低雷诺数通气空泡流进行了数值仿真，数值结果详细描述了与空气和水的剧烈混合相关的湍流运动，以及空泡流中不同特征结构之间的复杂相互作用。 
当气流中的局部静压低于蒸汽压力时，自然空化就会发生，从而导致气泡或空化的形成。空化效应可以大大降低表面摩擦阻力，使水下航行器能够实现高速行驶，特别是当气泡足够大并包围航行器在水中行驶的大部分表面（或整个表面）时。在实践中，由于当航行器从最初的非空化状态加速到空化状态时所面临的技术挑战，例如由抖振和表面变形振动引起的不稳定效应，自然空化效应很少产生大的空化。另一方面，人工通气（包括吹入不可冷凝气体以形成空泡）可以克服其中的许多困难。然而，深入了解通气空泡的内部结构及其与湍流的相互作用对于设计合适的高速水下航行器是必要的（Semenenko，2001）。
在通气空化过程中，空气通过后密封区（也称为封闭部分）从空泡中逸出。水和空气在这个区域发生剧烈的混合流动。以前的研究已经确定了两种主要的闭合模式，即回射流（Re-entrant Jet，RJ）和双涡（Twin Vortex，TV）。当重力效应不如对流效应显著时，就会出现RJ闭合模式。在这种情况下，空泡的后部充满了周期性排出的泡沫，这些泡沫呈甜甜圈状轴对称漩涡的形式（例如见Campbell & Hilborne 1958）。双涡管模式的特点是在空泡后部有一对涡流管，通常在重力效应占主导地位时出现（Cox & Clayden 1956）。在某些攻角条件下，通气空泡中的TVs具有与长椭球体中的TVs相似的几何形状。然而，在倾斜长椭球体情况之外的流动中，TVs是由涡流拉伸和倾斜过程产生的。而在空化流中，TVs是由重力效应产生的，这导致空泡在流向方向上上升，从而导致上下空泡表面之间的速度差以及随后的流动分离（参见Wang et al. 2021b）。大多数先前的研究侧重于在控制这些闭合行为的无量纲参数之间建立联系（例如：Epshtein 1973; Logvinovich 1973; Spurk 2002; Kinzel, Lindau & Kunz 2009; Wu et al. 2019）和了解导致不同流动条件下闭合模式变化的机制（例如：Semenenko 2001; Zhou et al. 2010; Kawakami & Arndt 2011; Ahn et al. 2017; Wu et al. 2019）。虽然在以前的研究中已经获得了大量关于闭合行为的知识，但很少有研究调查通气空泡中闭合部分附近的湍流动力学，这是本文的重点。闭合附近的湍流对于理解通气空化流中空气和水混合、泄气和涡脱落的机理非常重要。它也是空化流应用中湍流模型开发的物理基础。然而，尽管一些研究讨论了湍流，但大多数先前的工作集中在湍流和非定常行为的定性描述上，而没有定量分析空泡-湍流相互作用。特别是，很少有研究探讨涡流动力学的细节、湍动能的分布、能量转移收支或通气空泡闭合部分附近的雷诺应力分布。这些问题很难通过实验来研究，因为在空气-水混合区获得湍流动力学的非侵入式测量存在挑战。对于计算工作而言，准确模拟具有大密度比的多相流和大尺度范围的湍流运动是一项挑战。关于空泡-湍流相互作用的早期理论工作基于统计理论（如Plesset & Prosperetti 1977; Rood 1991）。后来的研究侧重于对空泡-湍流相互作用的更详细描述。例如，Callenaere等（2001）研究了RJ不稳定性，证明空泡密封区的压力梯度对封闭和局部湍流强度有很大影响。正如Arndt（2002）所述，空化会影响涡流动态并以复杂的方式产生湍流，而空泡是闭合区域附近的涡流源。最近，Gnanaskandan和Mahesh（2016）对圆柱体上的空化流进行了大涡模拟（LES）。他们发现空化可以减弱空泡尾流附近的湍流。Karathanassis等（2018）证实了这一结果。Koukouvinis等（2017）使用实验技术和数值方法分析了孔口内的空化流，发现孔口区域的湍流受到空化结构的影响。更具体地说，湍流在密集空化区域被抑制，而在空泡溃灭区域被增强。此外，他们发现，湍流在存在完善的空化时会受到抑制。Barbaca等（2019）研究了自由流雷诺数和弗劳德数对空泡闭合动力学的影响，并对通气流和自然空化流进行了比较。他们应用阴影成像和适当正交分解（proper orthogonal decomposition, POD）分析技术研究尾流，注意到自由流速度影响通气空泡流湍流混合区中脱落结构的大小和数量。Wang等（2021a）通过数值模拟研究了圆柱尾流中的通气空泡流。他们在湍流动力学方面的主要发现是，随着气体夹带量的增加，封闭部分附近的湍流强度降低，同时气泡数量减少。此外，先前的几项研究表明，相干流结构在通气空泡的湍流动力学中起着重要作用（例如：Dittakavi、Chunekar & Frankel 2010; Wang等 2018; Wu等 2021）。然而，我们仍然缺乏对通气空泡闭合部分附近的相干结构如何与湍流动力学相关联的了解，以及这些相干结构在泄气和涡脱落过程中所起的作用。为了填补这一知识空白，在目前的工作中，我们采用了一种谱本征正交分解（spectral proper orthogonal decomposition, SPOD）方法来识别通气空泡流中闭合区域附近的拟序结构，以通过模态分解方法有效地描述动力学。SPOD方法是最初由Lumley（1967，1970）提出的POD方法的变体，该方法提取一组本征模，以开发一组给定流量数据的最小二乘表示。Sirovich（1987）描述的POD形式之一已在文献中广泛使用（例如Aubry 1991; Berkooz, Holmes & Lumley 1993; Meyer, Pedersen & Özcan 2007; Hellström, Ganapathisubramani & Smits 2015）。这种POD方法旨在通过分解空间相关张量来识别空间正交模式，假设每个瞬时快照都是时间独立的实现，并被Towne，Schmidt & Colonius（2018）称为仅空间POD方法。相比之下，SPOD方法通过分解交叉谱密度张量来提取模式，每个模式都以特定的频率振荡。因此，SPOD方法解释了快照中的时间一致性（Towne等，2018）。SPOD方法已被用于识别各种流动环境中的相干结构。例如，Nidhan、Schmidt和Sarkar（2022）使用SPOD方法提取和分析圆盘湍流尾流中的相干结构，并揭示了它们与浮力、涡脱落和非定常内部重力波的关系。Abreu等（2020）利用SPOD方法分析了湍流管流中的相干结构。Li等（2021）使用SPOD方法研究了地面车辆附近的相干结构及其空气动力学特征。Ghate，Towne & Lele（2020）通过将SPOD方法与基于物理的富集算法相结合来重建湍流场，该算法使用时空局部化的Gabor模式来表示惯性范围内的湍流。Kadu等（2020）利用SPOD阐明了无约束同轴射流中的动态重要模式。Schmidt等（2018）使用SPOD方法复现了喷流中的湍流结构，并将不同位置观察到的波包与流动动力学联系起来。Nekkanti & Schmidt（2021）展示了SPOD方法在湍流射流环境中的几种应用，包括流场重建、去噪、频率-时间分析和预白化。在本文中，我们首次采用SPOD分析通气空泡流中的拟序流结构。我们的研究利用了直接数值模拟（DNS）方法，该方法显示了足够的分辨率来解析与闭合和夹带机制相关的关键动力学。随着计算能力的增强，计算流体动力学（CFD）已成为空化研究中获得三维流场详细描述的强大工具。除了上述Gnanaskandan & Mahesh（2016）和Wang等（2021a）的工作外，Uhlman & James（1989）和Abraham，James & Ivan（2003）使用边界元法对稳态空化流进行了数值模拟。Wang和Ostoja-Starzewski（2007）使用雷诺平均纳维尔-斯托克斯（RANS）模拟和大涡模拟研究了自然空化引起的空化流，在模拟中将两相流体视为相干流体系统。Senocak & Shyy（2002）使用两相RANS方法模拟了非定常空化流。Kunz等（2000）提出了一种预处理纳维尔-斯托克斯方法来模拟空化动力学。Lindau等（2015）基于有限体积法进行了RANS模拟和大涡模拟，以研究空泡中的脉动现象。kinzel（2008）通过结合level-set界面捕捉方法和重叠网格方法开发了一种可压缩空化流的数值方法，从而提高了模拟的精度和稳定性。这个数值工具也被用来研究空泡流的空气夹带过程和闭合模式。本文介绍了通气空化流的DNS研究，重点是空化和湍流之间的相互作用，包括两相湍流、内部结构和统计量的综合分析。我们通过分析通气空化流中的泄气、涡脱落和振荡过程来研究各种空泡结构中的湍流动力学。这项研究的目的是解决以下问题。（3）闭合区附近的特征结构在湍动能（TKE）输送中起什么作用？（4）究竟是什么流动过程控制封闭区域的泄气和涡脱落？本文的其余部分结构如下。首先，第2章介绍了控制方程、数值方法和模拟方法。第3章介绍了流场的时间平均统计数据（包括平均流量、雷诺应力和涡流动力学）。在第4章中讨论TKE及其在空泡中的收支。此外，第5章研究了泄气及其与相干涡流结构的相关性。在第6章中用SPOD方法分析了闭合区附近的相干结构。第7章讨论了低雷诺数DNS模拟的评估，以深入了解通气空化流。最后，在第8章中总结了结论。在这一节中，我们首先介绍2.1节中的问题设置和控制方程。然后，在2.2节中，我们描述了用于模拟空气和水流动并跟踪动态空气-水界面的数值方法。最后，我们在2.3节中给出了计算参数。01 问题设置和控制方程
在我们对通气空化及其湍流的机理研究中，我们考虑了一种计算设置，其中通气空化器位于矩形水洞内，从入口均匀进水。在当前模拟中，我们仅考虑主流速度不是很高的情况，因此无法实现自然空化状态。换句话说，在我们的计算设置中没有考虑自然空化效应对空泡的影响。图1（a）显示了我们坐标系中的计算域，x、y和z分别表示流向、垂直和横向方向。水以速度U∞沿x方向流过入口，通过空泡，并通过出口流出。在模拟中，在入口处施加狄利克雷边界条件，在出口处施加对流边界条件（Lodato，Domingo & Vervisch 2008）。
空化器左侧的蓝色箭头指示进水方向；绿色箭头指示空气流通方向 图1（b）显示了空泡示意图，显示了空化器和通气源。空化器位于xc=30dc处的中心轴上，是一个直径为dc、厚度为0.3dc的圆盘。我们使用直径为0.6dc的球形源区域，该区域位于圆盘后面1dc处，以在模拟中产生通气。这种数值处理方法克服了用有限的网格分辨率表示实际小空气通气配置的数值困难。更多的数值细节在2.2节中给出。在我们的实验和应用中，空气通过安装在支撑空化器的支柱上的排气口排出。通气口的形状和位置以及空气流通方向会影响空泡（Semenenko 2001）。因为目前工作的目的是阐明空泡的基本机制，所以复杂的物体和通气口设置可以用一个简化的模型来代替，而不会丧失这种机理研究的一般性。在模拟中，使用密度和粘度变化的流体系统同时处理空气和水。流体运动由不可压缩的纳维尔-斯托克斯方程和连续性方程控制，如下所示:在上述等式中，xi (i=1, 2, 3)表示笛卡尔坐标(x, y, z)，ui是流速的(u, v, w)分量，ρ是流体的密度，p是动压（此处作者应为笔误，应为去除重力项后的静压），gi (i=1, 2, 3)表示重力加速度的分量，σij=μ(∂ui/∂xj+∂uj/∂xi) (i,j=1, 2, 3)表示应力张量的分量，Ti (i=1, 2, 3)表示水气界面的表面张力的分量，μ是动力学粘性系数。02 数值方法
动量方程（2.1）使用二阶中心差分方案进行空间离散化。当我们采用这个方案时，注意到几个考虑因素。首先，这种方案易于实现，并且每一步所需的计算时间比其他方法少。此外，该方案与我们用于求解泊松方程的方法一起稳健地工作。此外，该方案的短模板允许我们实施浸入边界方法和耦合level-set和流体体积（CLSVOF）方法，并增加了数值鲁棒性。最后，该方案在数值上是无耗散的，这使得DNS能够在大范围的尺度上解析湍流运动。时间推进采用二阶龙格-库塔（RK2）方法。为了加强连续性方程（2.2），分步方法（Kim & Moin 1985）在RK2方法的每个子步骤中使用。压力通过求解泊松方程获得，如下所示:其中，Δt为时间步长，u*为投影法中的中间速度，但未在（2.1）中添加压力梯度。此外，是通气源中的发散常数，φ是Heaviside函数，在通气源内部等于1，否则等于0，ds=0.6dc是通气源的直径，u0是均匀空气通气速度（如图1a所示）。可以使用空气夹带系数Cq=Q/（U∞dc2）计算均匀空气通气速度u0，其中Q是空气通气率:方程（2.3）是一个具有大矩阵的泊松方程，因为它涉及三维网格点。如此大的矩阵的直接求逆在计算上是禁止的。因此，这个矩阵需要用迭代法求解。此外，模拟必须在大型并行计算机上进行，这增加了计算成本。可移植的可扩展科学计算工具包（PETSc）包已在许多先前的研究中得到应用（例如，参见Balay等2017; Yang, Deng & Shen 2018; Gao, Deane & Shen 2021b），并已证明在求解泊松方程时计算效率高且稳健。因此，本研究采用了该软件包。空气-水界面使用CLSVOF方法获取（Sussman & Puckett, 2000）。level-set函数ψ描述了到空气-水界面的距离，在水和空气中分别为正值和负值。ψ的零值对应于界面，ψ的梯度的范数是一致的，即其可以通过应用二阶本质非振荡方案（Shu & Osher 1989）进行空间离散化，并通过二阶算子分裂方法进行时间积分。采用重新初始化程序（Min 2010）来确保对流ψ是（2.6）所述的带符号距离函数。气相和水相通过其密度和粘度来识别，密度和粘度根据level-set函数ψ确定其中H(ψ)是Heaviside函数，ρw和ρa是水和空气的密度，μw和μa分别是水和空气的动力学粘度。在模拟中，这些值在空气-水界面附近被平滑为这里的ϵ是扩散区的一半厚度。根据我们之前的研究（Yang et al. 2018），设置为最小网格大小的三倍。根据颜色函数φ捕获界面的流体体积（VOF）法用于强制每个流体相的质量守恒。这里，φ是网格单元中水的体积分数，在水中等于1，在空气中等于0。φ的控制方程为这是通过使用守恒性的VOF方法解决的（Weymouth & Yue 2010）。Sussman，Smereka & OSher（1994）提出的CLSVOF方法用于耦合ψ和φ。在每个时间积分子步骤中求解对流方程（2.7）和（2.11）之前，level-set函数ψ用于确定每个网格单元中界面的法线方向。法线方向的系数通过最小化从周围节点到界面的距离的加权和来计算。在获得每个网格单元中的界面函数后，使用分段线性界面计算算法计算φ的通量，并将其用于φ的对流。然后φ用于调整界面附近对流ψ的值，以保持质量守恒。有关CLSVOF算法及其在各种多流体问题中的应用的详细信息可在我们小组以前的论文中找到，包括Hu等（2012）、Yang等（2017、2018）和Gao等（2021a）。为了捕捉水下空化器对气流的影响，采用了浸没边界法（Mittal & Iaccarino 2005; Kim & Choi, 2019）在流动求解器中使用（Cui等 2018），其中在空化器表面附近的网格点上施加适当的体力，以加强无滑移边界条件。浸没边界法可用于在笛卡尔网格中而不是在贴体网格中解析固体，从而简化算法并减少计算时间。图2显示了我们模拟中的典型瞬时空泡。图中显示了几个通气空化组件，包括界面（空气-水边界）、前部（空泡主体）和封闭部分（空泡后部，空气从中逸出，两相在此剧烈混合）。为了进一步说明这些动态，提供了一部补充影像。它提供了实时的视觉呈现，特别突出了剧烈的混合过程和空气包从空泡的封闭部分逸出（补充影像请参见https://doi.org/10.1017/jfm.2023.431）。
03 计算参数
50dc×50dc×50dc（长×高×宽）计算域使用拉伸笛卡尔网格进行离散化。根据可用的计算资源，计算域的长度被扩大到最大可能范围，以确保出口边界条件对通气空化的影响不明显。我们使用不同分辨率的网格进行了三次模拟，分别表示为M1、M2和M3。M1、M2和M3的网格点分别为400×96×96、800×192×192和1600×384×384。网格在空化器周围细化，并沿(x-xc)/dc<4和(x-xc)/dc>18方向逐渐变粗（如图3）。表1根据最小单元尺寸Δmin、时间步长Δtf、网格数量、用于并行计算的内核数量以及用于模拟的总内核小时数列出了网格分辨率。
其中红色虚线矩形表示细化网格区域，绿色矩形表示圆盘空化器 出于说明目的，每个方向上每四个网格点中只有一个被绘制 M1模拟直到空泡的尺寸在大约tc=700dc/U∞处统计稳定。继续模拟tc=500dc/U∞以获得流量统计数据。然后，通过将t=tc时的M1模拟流场数据插值到具有更高分辨率的网格中来初始化M2和M3模拟。M2和M3模拟的流量统计数据分别在瞬态期t=3.8dc/U∞和t=5.2dc/U∞后进行计算，在此期间，流量根据新的网格分辨率调整到统计稳定状态。M2和M3模拟的流量统计数据分别在t=480dc/U∞和t=440dc/U∞的持续时间内进行计算。
时间步长根据稳定性标准设定（Ling等 2019）其中c是Courant–Friedrichs–Lewy（CFL）数，Δtconv是对流项的时间步长限制，tvis是粘性项的时间步长限制。以模拟M3为例。CFL数c设置为0.4，Δmin/dc=0.01，ΔtconvU∞/dc=0.005，ΔtvisU∞/dc=0.01。因此，我们将时间步长设置为Δt=2×10-3，以确保计算稳定性。通过与文献中的数据和网格收敛测试进行比较，验证了数值方法和网格分辨率选择方法，详情见附录A。通气空化流由几个无量纲参数表征，包括:空化数σc=(p∞-pc)/(0.5ρwU∞2)，其中p∞和pc分别为来流压力和空泡内部压力；通气系数Cq=Q/(U∞dc2）；弗劳德数Fr=U∞/√gdc，其中g为重力加速度；雷诺数Re=ρwU∞dc/μw；阻塞比B=dc/DT，其中DT代表水洞横截面的水力直径；以及韦伯数We=ρwU∞2dc/σT，其中σT为表面张力系数。请注意，σ的定义仅考虑了通气效应，这与Ahn等（2017）的先前研究一致，因为与通气气体相比，自然空化产生的水蒸气可以忽略不计。在这项研究中，我们感兴趣的是重力效应可以忽略不计的场景；因此我们将Fr设为∞。在这种情况下，流场是轴对称的，闭合模式是RJ。因为大部分空泡的曲率很小，表面张力的影响可以忽略不计，所以We设置为∞。其他无量纲参数设定为Cq=0.14、Re=1000和B=1.8%。选择相对较小的雷诺数，在不引入湍流模型的情况下通过DNS解析空泡中的所有湍流涡流。尽管在实践中雷诺数要高得多(105-107)，但目前的DNS方法仍然是一种有用的研究工具，可以深入了解流动物理的基本动力学（Moin & Mahesh 1998）。在之前的研究中，我们发现B=1.8 %足够小，可以防止阻塞效应。密度比设置为ρa/ρw=0.0012，动力学粘度比设置为μa/μw=0.0154，这是空气和水的实际值。尽管模拟是在笛卡尔网格上进行的，但必须将数据从笛卡尔坐标（即x，y，z）转换到圆柱坐标（即r，θ，x）来分析流场，因为Fr=0（此处作者笔误，应为∞），所以流场在统计上是轴对称的，如2.3中所述。转换按如下方式执行:圆柱坐标系中的矢量分量，即（fr, fθ, fx），可以根据笛卡尔坐标系中的分量，即（f1, f2, f3）计算如下:在模拟达到完全湍流和统计稳定状态后，变量f的平均量（表示为f（平均））可以通过在方位角方向和时间上取平均值来确定为其中i、j和k是x、y和z方向上的离散网格点的索引，以及在上面的方程式中，c是一个被设置为1.9的常量，以确保对于所有考虑的r值，都有： 
并且f不被过度平滑。统计平均化是在持续时间80dc/U∞上执行的，因为在此期间流体颗粒可以通过计算域。为了考虑水相和气相之间的密度差异，采用密度加权平均(也称为Favre平均)来分析结果。Favre平均值表示为因此，波动被定义为
01 平均流动
图4显示了叠加在平均流动流线上的水的平均体积分数的等高线，这些等高线是根据时间和空间平均(3.6)获得的。流线表明在空泡内存在一个回流区，在该回流区内，空气在气液界面附近向前运动，但在空泡中心轴附近向后运动。当φ（平均）=0.5时，空泡的平均轮廓由平均体积分数的等值线定义。图5显示了空泡闭合时水的平均体积分数和相应的等值线的等高线。封闭区域内的φ值在0.2~0.8范围内，表明该区域内空气和水的混合很剧烈。此外，由于外部水的射流，轮廓线在中轴线附近呈凹形。这一现象的基本物理学将在第3章的其余部分以及第4章和第5章中进一步讨论。
图4：水的平均体积分数φ（平均）的等高线和平均流场u（平均）的流线
图5：水的平均体积分数φ的等高线和空泡闭合区的相应等值线
空气相的流向速度和径向速度定义为（Kinzel等 2009）：
其中，ur是径向速度，根据(3.3)计算。它们描述了空泡内空气的流向和横向运动，因为如果位置在气相，它们等于局部速度，在水相等于零，在有两个相的网格单元中等于(1-φ)u和(1-φ)ur。
（a）Ua(平均)等值线，（b）U a(平均)沿中轴分布，（c）Ur,a(平均)等值线，（d）区域4.5<(x-xc)/dc<8.0和0<r/dc<1.4的Ua(平均)等值线，以及（e）Ua(平均)在区域11.5<(x-xc)/dc<14.0和0<r/dc<1的等值线。(a, c)中的黑色矩形和绿色半球分别表示空化器和通气源。黑色虚线对应于φ(平均值)=0.5的等值线。(a, c)中的垂直实线表示图7中绘制横截面轮廓的位置。绿色虚线矩形对应于(d, e)的位置。 图6(a,b)显示了腔内的ua(平均)和ur,a(平均)。图6(a)表示腔体内部的红色 区域和蓝色 区域分别发生向前和向后运动。图6(b)显示，与流向速度相比，气流的径向速度可以忽略不计，除了靠近空泡头部和后部的区域。红色 区域的φ(平均)=0.5等高线表示空气向外移动并接近空泡头部的界面，然后向内向靠近空泡后部的中轴线移动。这些结果与图4所示的流型是一致的。图6(c,d)显示了腔的中间和闭合区域中的ua(平均)的轮廓。图6(c)显示两种流动结构在空泡的中部占主导地位：剪切层(SL)，对应于界面附近的红色等高线区域，以及回流区(RA)，对应于流向中轴附近的蓝色 区域。在SL中，流动以高局部剪切为特征。在RA中，流体向后向空化器方向移动并减速。如图6(e)所示，在闭合区域附近的空泡外可以观察到另一种特征流动结构，称为喷射层(JL)。喷流的一支向后运动，以两相混合流体的形式重新进入气泡，而另一支喷流向前运动进入尾流。SL和RA的结构在以前的研究中已经讨论过(Callenaere等 2001；Spuk 2002；Kinzel等 2009)。然而，据我们所知，以前还没有研究过通气空化流动中的JL。此外，SL、RA和JL与空泡内的湍流之间的关系尚未阐明。在本研究中，这些结构在湍流产生和传输中的作用在第4章中进行了研究。此外，在第5章中，讨论了它们与泄气和涡脱落过程的关系。
(a)流向速度ux,a(平均)；(b)径向速度ur,a(平均) 图7(a,b)显示了在图6中垂直线标记的流向位置处流向和径向分量的平均气相速度分布。图7(a)表示最大流向再循环速度，即最大负ua(平均)，出现在空泡的中心轴处。此外，其大小沿流向方向增大。正ua(平均)的峰值出现在空泡的气液界面附近。该值增加，直到(x-xc)/dc∼6处的最大空泡半径，然后减小。图7(b)显示，平均径向速度ur,a(平均)的峰值在该点之前减小，然后以相反的符号增加。此外，图7显示空气相流向运动主导腔内的径向分量。在文献中，气相流向速度Ux,a也被称为空气夹带的局部速率(Kinzel等 2009)。
图8显示了平均空化数σc(平均)沿中心轴的分布。在封闭区域附近，由于JL区域内部流动停滞导致的压力增加，σc(平均)突然减小。然而，σc(平均)在腔内保持不变。02 涡动力学
由于难以获得空泡内部的测量结果，因此仅进行了少数关于空泡涡动力学的实验研究（例如Yoon等 2020）。然而，涡动力学在泄气、通气和湍流产生中很重要。本文研究了泡内各部分的特征涡结构和相应的涡动力学，以阐明空泡的物理机制。
图9显示了计算域中通过中面的瞬时空间涡度ωθ的等高线。上方的蓝色 区域和下方的红色 区域指示了空泡前缘的两个SL结构，涡度主要是由水和空气相之间的切变产生的。在封闭区，这两种SL结构相互作用，产生较小的涡旋结构。此外，这些小漩涡中的一部分离开空泡进入尾流。在再循环区域，即空泡的前部，可以观察到红色和蓝色交替的涡旋结构。Yoon等 (2020)也在实验上观察到了空化器附近的涡。然而，这些涡是由一个不同的过程产生的，即射流和回流之间的相互作用。机制的不同是由于本模拟和Yoon等 (2020)的实验之间的通气设置的差异引起的。Yoon等 (2020)将射流的方向设置为反向流动，这会导致两种流动之间的相互作用，从而产生漩涡。相反，在目前的模拟中，我们使用球形源将空气通气到空泡中，这不会导致通气空气和反向流动之间的强烈相互作用。模拟中观察到的涡旋主要产生在闭合区，并由再循环流向后对流。涡度的平方ζ=(ωx2+ωy2+ωz2)/2表示单位质量的涡度。在这项工作中，这个量被用来评估涡旋结构的强度，以研究腔内的涡动力学。其中ζx=ωx2/2是由横截面上的螺旋运动引起的顺流拟能，ζyz=ωyz2=(ωy2+ωz2)/2表示横向拟能且主要与x-r平面上的剪切运动有关
图10：平均拟能及其分量的云图 (a)ζ(平均); (b)ζx(平均); (c)ζyz(平均)
虚线表示φ(平均)=0.5时的等值线; 垂直实线表示图11中横截面轮廓的位置
图10显示了平均拟能ζ(平均)及其两个分量ζx(平均)和ζyz(平均)的云图。图10(a)表明，在空泡的前部和封闭区的SL中发生了强烈的涡旋运动。图10(b)显示强流向涡度只出现在封闭区。此外，图10(c)表明，强烈的横向涡旋不仅出现在空泡的闭合区，而且也出现在前部的SL中。在闭合区域，横向拟能大于流向拟能。空泡内横向涡动的主要原因是气液界面的强剪切作用。封闭区域内强烈的流向涡动表明，由于该区域内两种流体的强烈混合，该区域内的流动转变为更加三维和各向同性的状态。图11：不同下游位置的平均拟能其分量的分布：(a)ζ(平均)、(b)ζx(平均)和(c)ζyz(平均)
为了阐明空泡中涡旋动力学的基本物理机制，我们研究了沿径向的平均拟能分布。图11显示了泡内不同流向位置(参见图10中的垂直线)的平均拟能(ζ)及其两个分量(ζx(平均)和ζyz(平均))。横向拟能ζyz(平均)的峰值出现在空泡界面附近，剪切效应较强。峰值ζyz(平均)在前面大致保持不变(比较(x-xc)/dc=2、4、6、8、10)，并在闭合区域((x-xc)/dc=12)变大。然而，流向拟能的峰值沿下游方向增大，直至闭合区域。产生上述变化的原因是横向涡动是由两相切变效应引起的，两相切变效应使前半部分的SL中的强度保持不变，而流向涡动则归因于湍流的发展，湍流沿流向逐渐增强。为了研究涡旋动力学和确定空泡-湍流相互作用中的调节机制，我们研究了拟能收支方程，该方程表示为这个方程是从变密度的Navier-Stokes方程(2.1)导出的。涡旋动力学由三项控制，即拉伸效应项Ws、斜压效应项Wb和粘性效应项Wv在(3.16)中。这些效应如图12所示。伸展和斜压效应是产生拟能的原因，而粘性项则使拟能消散。图12显示，这三种效应在空泡的前部都很弱，在闭合区很强。
图12：拟能方程中平均收支项的云图(3.16)：(a)伸展效应项Ws，(b)斜压效应项Wb，和(c)粘性效应项Wv。虚线表示φ(平均)=0.5时的等值线;垂直实线表示图13中横截面轮廓的位置。
图13：收缩方程（3.16）中不同流向位置的平均收支项的剖面：(a)拉伸效应项Ws(平均)，(b)斜压效应项Wb(平均)，以及(c)粘性效应项 Wv(平均) 此外，我们还得到了三个收支项的径向分布。图13显示了这三项在不同下游位置的云图。斜压项和粘性项的大小是相当的，并主导着伸展项。斜压效应在涡旋的产生中起着至关重要的作用，而斜压项和粘性项几乎相互抵消。此外，这三项的峰值都出现在空泡内的气水界面附近，并且由于气水两相相互作用的增强，沿下游方向有增大的趋势。03 雷诺应力
雷诺应力在动量输运以及湍流脉动与平均流的相互作用中起着重要作用，而平均流是TKE产生的主要原因。以前的研究已经研究了片空化和附着式空化过程中的雷诺应力分布(Gopalan & Katz 2000；Iyer & Ceccio 2002)。然而，对主动通气空化过程中的雷诺应力分布知之甚少。目前的模拟数据可以用来解决这一知识差距。在空气-水两相流动中，雷诺应力是以可变密度来定义的。变密度流的Favre平均动量方程写为是雷诺应力张量。在RANS模拟中，雷诺应力需要根据平均量闭合。为了开发描述未来空化流的模型，我们研究了雷诺应力的结构特征和空间分布。由于所研究的空化流的轴对称性，雷诺应力分量τRxθ和τRrθ为零。因此，在下面的讨论中只考虑雷诺应力张量的四个独立分量，即τRxx、τRrr、τRθθ和τRxr。
图14：雷诺应力分量的等值线：（a）τRxx、（b）τRrr、（c）τRθθ和（d）τRxr。黑色虚线表示φ(平均值)=0.5的等值线。垂直线表示图15中绘制横截面轮廓的位置。
图14显示了雷诺应力张量的四个分量的云图。为了更清楚地可视化分布，我们还在图15中描绘了几个下游位置的雷诺应力分量剖面。如图所示，在空泡的前部，τRxx是雷诺应力的主要分量，而其他分量则可以忽略不计。图15（a）表明τRxx的最大值出现在空气-水界面处，两相之间发生强相互作用，峰值在下游方向增加。图15（b，c）显示，τRrr和τRθθ在空泡的闭合区域相似，表明该区域的流动变得更加各向同性。此外，图14（d）显示τRxr大多为正。附录B中的象限分析表明，正τRxr分量与雷诺剪切应力Q1和Q3象限Q2和Q4象限的主要贡献有关。同样，τRxr的峰值往往出现在两相界面附近，并沿流向增加。
图15：不同流向位置的雷诺应力分量剖面：（a）τRxx、（b）τRrr、（c）τRθθ和（d）τRxr 综上所述，根据雷诺应力分析可以得出以下结论。首先，雷诺应力的流向分量的峰值出现在空气-水界面附近，并且该分量在大部分空泡中占主导地位。其次，雷诺应力在闭合区附近变得更强，并且对角线分量具有可比性，表明该区域的流动变得更加湍流和各向同性。TKE动力学在湍流建模中尤为重要。在本多相流研究中，基于Favre平均值以及（3.11）和（3.13）中定义的波动，TKE定义为
图16：Favre平均TKE的等值线，q(平均)。绿色实线是q(平均)的等值线，
黑色虚线对应于φ(平均)=0.5的等值线
湍流主要发生在闭合区，往往集中在JL结构中。图16绘制了该区域Favre平均TKE q(平均)的等值线。TKE等值线在定性上类似于一组椭圆。它们的中心位于JL的停滞点附近，TKE在那里达到峰值。强烈的湍流是由该地区空气和水的剧烈混合引起的。为了进一步评估腔体中湍流的各向同性，我们考虑了TKE分数函数（Liu & Xiao 2016），其表示为该分数是与流向波动相关的TKE与与横向波动相关的TKE之比的两倍。如果fq∼1，则湍流达到近似各向同性状态。图17显示了与图16所示的同一闭合区域中fq的轮廓。峰值出现在JL的剪切层和停滞点附近。峰值大于3.0，表明在SL和JL中，TKE的流向分量都比TKE的横向分量强得多。
图17：TKE分数函数的等值线，fq。绿色实线是fq的等值线，黑色虚线对应于φ(平均)=0.5的等值线 为了阐明TKE的产生和运输机制，我们基于密度加权平均研究了其收支方程，如下所示：其中Πc是平均流动对流输运项，Πd是产生项，Πtt是湍流输运项，Πp是压力输运项，Πvt是粘性输运项，ϵ是粘性耗散项。数值结果表明，Πtt和Πvt可以忽略不计（此处未绘制）。因此，在下面的讨论中，我们将重点放在其余项上。
图18：粘性耗散项云图
绿色实线是ϵ的等值线，黑色虚线对应于φ(平均)=0.5的等值线
图18显示了耗散项的云图，其中闭合区域中TKE居于最大值（图16）。如图16所示，从定性上，云图分布与TKE相似，特别是在JL区域。最大耗散发生在停滞点附近。众所周知，它总是正的，并将TKE转化为内能。
图19：对流输运项Πc云图
箭头绿线表示流线，黑色虚线对应φ(平均)=0.5等值线，黑色实线表示TKE等值线
图19显示了对流输运分量Πc的云图。SL中的蓝色 区域表示能量通过平均流从该区域输出，而JL区域中的红色 区域表示平均流将能量输运到该区域。这种TKE输运现象可以解释如下。我们可以将Πc改写为Πc=-ρ(平均)u(平均)k∂q/∂xk=−v·∇q。此外，图19显示了TKE的云图。在这里，∇q垂直于这些等值线并指向停滞点。流线被绘制为带有箭头的绿线。在SL中，平均流量向停滞点移动，导致负Πc(−v·∇q<0)。因此，TKE输运出SL区域。然而，在JL中，平均流量远离停滞点，导致正Πc(−v·∇q>0)。因此，TKE被对流到两个JL分支中。图20（a）绘制了TKE生成项Πd的云图，在SL和JL区域显示出不同的特征。红色 区域的正等值线值表示SL中的Πd将动能从平均流传递到湍流，而蓝色 区域的负等值线值表示能量向相反方向传递。为了阐明这种能量转移过程背后的物理机制，我们将Πd分解为三个部分： 
其中Πdx表示流向拉伸效应产生的产量，Πdt是横向运动的产量，Πds是剪切效应产生的产量。图20（b）显示，SL和JL中的Πdx均为负值。图20（c）表明Πdt为正，并集中在JL中。图20（d）显示Πds主要为正，并集中在SL中，请注意，Πd在SL中为正，在JL中为负。因此，Πds在SL中占主导地位，而Πdx在JL中占主导地位。为了进一步评估上述项的相对重要性，我们计算了不同项的大小的分数函数，如下所示
图20：TKE生产项及其组成部分的轮廓：（a）Πd，（b）Πdx，（c）Πdt和（d）Πds。
黑色虚线表示φ(平均)=0.5的等值线
在SL和JL中，与fdx和fds相比，fdt可以忽略不计（未绘制）。fdx和fds的轮廓如图21（a，b）所示：fdx集中在JL中，而fds在SL中占主导地位。这些结果可以通过腔体闭合区域的雷诺应力和平均应变率的分布来解释。在SL中，剪切效应导致正而强的∂u(平均)x/∂r。此外，雷诺应力τRxr是正的且不可忽略的。因此，TKE的生产率以正Πds分量为主。然而，在两个JL分支中，射流增加了拉伸应变率分量，即射流产生负∂u(平均)r/∂r和正∂u(平均)x/∂x分量，幅度很大。请注意，流向雷诺应力在JL中占主导地位（见图17）。因此，由于τRxx为负，而∂u(平均)x/∂x为正，因此JL中的TKE生产率由负Πdx分量主导。
图21：Πd分量的分数函数的轮廓：（a）fdx和（b）fds。
黑色虚线表示φ(平均)=0.5的等值线。
其中Πpm是平均压力梯度引起的输运，Πpt是压力-速度波动相关性引起的输运，Πpd是膨胀效应引起的输运。请注意，由于两相的不可压缩性，Πpd会降低。Πpm和Πpt的云图如图22（a，b）所示。图中显示，SL和JL的Πpm均为正值，增强了腔体闭合区的TKE。另一方面，Πpt在SL中基本上是负的，并且倾向于将TKE从SL区域转移到JL区域。
图22：压力传输项的云图：（a）Πpm和（b）Πpt。
黑色虚线表示φ(平均)=0.5的等值线
图23：平均流量、SL湍流和JL湍流之间的能量传输示意图 总结第4章，图23说明了SL和JL区域湍流与平均流量之间的能量转移过程。上述结果表明，SL和JL结构在TKE生成和运输过程中具有明显的特征。在SL中，由于剪切效应项Πds，TKE因生产项而增加，而在JL结构中，由于流向拉伸效应项Πdx，TKE减小。压力传递项Πpm倾向于增强JL和SL中的TKE。此外，平均流量对流输运项Πc和压力输运项Πpt有利于将TKE从SL输送到JL。在通气空化中，一些空气通过闭合区域从空泡中逸出，水和空气之间发生剧烈的流动混合。Wu等（2019）进行了实验，并使用PIV方法定量研究了RJ和TV闭合模式下发生的气体泄漏过程。此外，他们还讨论了RJ和TV模式转换期间气体泄漏的变化。他们的结果表明，气体泄漏是由于内部和外部流动之间的相互作用而发生的。根据目前的模拟结果，我们进行了多次观察，以阐明腔体RJ闭合区的流动特性和动力学。首先，仿真数据表明，封闭区域的泄气过程与涡脱落过程具有一一对应关系，如流动的瞬时可视化所示；以前的研究（Shao 等，2017; Karn，Arndt和Hong 2016）已经表明，RJ模式通过环形涡旋泄泄气体。然而，这些工作并没有将气体泄漏与涡脱落过程联系起来。图24（a，b）绘制了VOF函数φ=0.5和应变速率张量Q=2的第二个不变量的等值面，前者根据压力的大小着色，后者根据φ的大小着色。泄漏的空气-水界面表现出与Q（Haller 2005）特征的涡流结构相同的模式。此外，只有少数先前的实验室研究（Karn 等，2016 ;Wu 等，2019）确定了气泡离开空泡的位置以及相关的涡脱落。从我们的仿真结果来看，流动可视化表明，在闭合区域的JL结构内发生了泄气和涡脱落。图24（a）中的颜色说明了空泡表面的压力分布。泄漏的气流结构与腔体的断开和破裂发生在高压区域附近，这与JL结构内部的停滞区域相对应。图24（b）中的涡旋结构根据等值面上的φ值进行着色。结果表明，涡旋管和涡环发生在剧烈的空气-水混合区域，即JL区域。
图24：空泡封闭区域泄气和涡脱落的瞬时结构
（a）φ=0.5的等值面，根据压力的大小着色。（b）应变率张量（Q=2）的第二个不变量的等值面，根据φ的大小着色。
此外，我们观察到腔体闭合区的泄气和涡脱落过程不是随机的，并且这些过程表现出具有内在特征频率的半周期模式。由于测量限制，以前的实验工作很少研究闭合模式的光谱和振荡模式。Jiang，Shao和Hong（2018）对空泡长度进行了频谱分析，但没有发现主要振荡频率。在本研究中，我们的三维仿真结果说明了闭合模式下不同结构内振荡的更多细节。
图25：径向速度u和VOF函数沿中心轴φ的频谱：（a）log10( E(Ur) )和（b）log10(E (φ) ) 为了研究涡脱落与泄气过程之间的关系，我们采用时间傅里叶变换得到了径向速度E（ur）和VOF函数E（Φ）的频谱。沿中心轴的10<（x−xc）/dc<15的结果如图25所示。图25（a）和图25（b）都显示了一个重要的特征频率，即St≈0.1，其中光谱达到其局部峰值，表明泄气和涡脱落过程具有相同的特征频率。此外，光谱更集中在12.5<（x−xc）/dc<14区域，其中空气和水之间的湍流混合剧烈，导致该区域内的强烈湍流。泄气频率也可以通过将总气体体积定义为固定控制容积（CV）中的时间函数来量化，该控制容积以12.5<（x−xc）/dc<14和r/dc<0.8为界。CV中空气总体积Va（t）表示为
图26：区域12.5<（x−xc）/dc<14，r/dc<0.8区域气泡体积的归一化自相关函数 图26显示了RV(Δt)与时间间隔Δt的关系，由dc和U∞归一化。峰值出现在Δt=11.1dc/U∞，对应于Strouhal数St=0.09。该值接近径向速度和VOF函数的频谱中的峰值（图25）。换言之，泄气和涡脱落过程以相同的频率发生，并且这两个过程在通气空化流的封闭区域具有很强的相关性。总之，图24-26展示了有关泄气和涡脱落过程的几个特性。首先，涡流结构的产生与气泡的泄漏有关。这两个过程都发生在停滞的高压区域，空气和水相之间有剧烈的湍流混合。其次，闭合区的涡旋脱落呈现半周期性模式，在能谱中可以观察到内在特征频率。第三，泄气和涡脱落过程具有相同的特性频率。基于上述结果，我们总结了图27草图中腔体闭合区域的泄气和涡脱落过程。其潜在机制可以描述如下。空气通过SL向前对流，然后终止，这伴随着SL中涡度的产生，如图中的第1部分所示。由于腔体前部的持续通气，空气不能停留在SL区域的末端，因此进入封闭区域的JL，空气和水完全混合，如图27中的第2部分所示。最后，泄气并离开空泡，产生复杂的涡流结构，这些涡流结构向下游对流，如第3部分所述。
第3部分：泄气伴随着复杂大规模涡流结构的形成。
在本节中，我们采用SPOD方法（Towne 等 2018）来提取和分析通气空泡流闭合区域附近的相干结构。为了执行SPOD分析，将以恒定时间间隔获取的流场q的Ns快照分解为沿时间方向的一组块，如Q=[q(1), q(2), ..., q(i),..., q(Nfreq)]，其中i是块的索引，q=[ux, ur, uθ]是圆柱坐标中速度分量的向量，Nfreq是每个块中的快照数量（Nidhan 等 2022）。为了改进我们的光谱分析，我们在Schmidt & Colonius (2020)的处理之后，在相邻块之间包括Novl重叠快照。因此，此分解后的块总数可以计算为然后，在每个模块中执行时间傅里叶变换，以获得频域中的投影，如下所示其中，i是虚部单位，ω=2πSt是角频率，根据Strouhal数St=dcf/U∞计算。频率为ω的Nblk块的集 合表示为Qω=[qˆ(1)ω, qˆ(2)ω, ..., qˆ(Nfreq)ω]。傅里叶变换后，我们通过特征值分解计算出对应于频率ω的SPOD特征值和特征向量，如下所示其中W是加权矩阵，Λω=diag( λ(1)ω, λ(2)ω, ..., λ(Nblk)ω)是一个对角矩阵，其中包含根据从j=1到Nblk的能量递减排列的特征值。第j个SPOD模式Φ(j)ω计算公式为：采用基于能量的内积，通过将权重矩阵W定义为梯形积分矩阵来计算特征值和特征向量。因此，范数可以计算为其中D是积分域。有关SPOD配方的更多详细信息，请参见Towne等（2018）和Nidhan等（2022）。对于本工作的SPOD分析，我们使用N=4096个平面的连续快照，在z=0时跨越9.5<(x−xc)/dc<15.5和0<(y−0.5Ly）/dc<1.3，其中包括通气腔的封闭区域。采样快照的时间间隔为Δtdc/U∞=0.1。我们选择频率数为Nfreq=512，并在模块之间实现50%的重叠，即Novl=256，以平衡改进的频谱收敛和减少傅里叶变换中每个窗口之间的频谱纠缠。在傅里叶变换之前，每个数据块都乘以汉明窗口函数，以有效减少由于频谱泄漏引起的误差（Schmidt & Colonius 2020）。为了研究不同位置的特征谱特征的变化，我们在四个连续的域（D1、D2、D3和D4）中进行了SPOD分析，如图28所示。
图28：水φ(平均)的平均体积分数的等值线和闭合区附近平均流场的流线
图中显示了用于SPOD分析的四个子结构域D1、D2、D3和D4以及整个结构域Dt。
黑色虚线是φ(平均)=0.5的等值线
图29：不同域中能量的累积分数，Cj
（a）作为模态指数n的函数;（b）作为St的函数。图显示了四个考虑域，即D1、D2、D3和D4
我们通过研究SPOD模式对涨落能量的相对贡献来开始分析它们。图29（a）显示了能量的累积分数Cj(n)，作为四个下游域D1、D2、D3和D4中SPOD模态指数n的函数，计算公式为：不同子域间曲线的比较表明，领先的SPOD模式贡献的相对能量随着下行距离的增加而增加。此外，我们观察到D1和D2的曲线基本上坍缩到一条曲线上，而D2和D3的曲线差异很大，表明D2和D3之间发生了过渡。这一结果与我们在图25中的观察结果一致，即D3包含空气和水剧烈混合且湍流强度强的区域。此外，图29（a）表明，第一种SPOD模式贡献了总能量的很大一部分，分别贡献了D1和D2中约60%的能量，以及D3和D4中约70%和80%的能量。这些结果表明了低秩行为，这些行为用于描述与领先模式相关的物理机制占主导地位的场景（Schmidt 等 2018）。SPOD分析表明，在各种流动动力学中可以观察到低秩行为，包括湍流分层尾流（Nidhan 等 2022）、喷流湍流（Schmidt 等 2018）和湍流预混火焰（Brouzet 等 2020）。在本研究中，低秩行为在通气空化流的封闭区附近普遍存在，图29（a）中的结果表明，随着封闭区下游距离的增加，低秩行为变得更强。
SPOD模式贡献的能量累积分数作为无量纲频率的函数，由St测量，计算公式为如图29（b）所示。图中显示，在下游方向，0<St<0.4范围内的Cj(St)先增加（D1到D3），然后减少（D3到D4）。然而，较高频率范围（0.4<St）的Cj(St)单调增加。该结果表明，如第5节所述，闭合区频域中的能量重新分布是由SL和JL结构之间的相互作用引起的。此外，图29(b)显示曲线在特征频率St=0.1附近具有拐点，这在D3和D4中更为明显。该结果表明，围绕特征频率振荡的SPOD模式贡献了大量的能量。
图30：在四个域中获得的SPOD特征谱：（a）D1，（b）D2，（c）D3和（d）D4。深色到浅色对应于λ(i)中递增的模型指数i。λ(1)和λ(2)之间的红色阴影表示它们的差异。
接下来，我们讨论SPOD特征谱。图30显示了不同下游域中的SPOD特征谱λ(i)。图中的红色阴影表示每个域中第一和第二模式之间的显着差异。这一结果与我们之前对低等级行为的观察一致。此外，低频范围内的差异往往更大，表明低频下的低秩行为更为明显。此外，该图显示了大约St=0.1处的明显光谱峰值，也可以在我们的其他光谱结果中观察到（图25和26）。该特性频率与该频率下相干结构的振荡动力学相关。此外，通过比较图30（c，d）和图30（a，b），我们发现D3和D4中不同的光谱峰高于D1和D2中的光谱峰。出现这种差异是因为D3和D4的湍流强度强于D1和D2。基于上述分析，我们从图29和图30中总结如下。（i）主要SPOD模式的相对贡献随着下游距离的增加而增加。（ii）SPOD模式显示低秩行为。第一种模式占主导地位，在低频范围内更为明显。（iii）SPOD特征谱中存在明显的谱峰值。
图31：径向速度Φ(1)r(x, r, St)/||Φ(1)r(x, r, St)||∞的归一化第一SPOD模态的实部等值线，在域Dt中：（a）St=0.10，（b）St=0.20，（c）St=0.31，和（d）St=0.92。黑色虚线是φ(平均)=0.5的等值线。（c）中的黑色实线表示观测到的波包，包括SLKH、JLr、JLout和JLin。
为了分析闭合区附近不同频率振荡的主要相干结构，我们研究了径向速度Φ(1)r(x, r, St)的第一个SPOD模式的实部，该模态是通过在整个域Dt中进行SPOD分析获得的，并绘制在图28中。图31中的结果与通过相位匹配在四个子域中获得的模式获得的结果一致，并提供了全局视图。根据第3章中的分析，在空泡闭合附近可以观察到两个特征性结构，即剪切效应占主导地位的SL和由强压力梯度引起的拉伸效应占主导地位的JL。根据图31，我们可以识别位于SL和JL附近的更精细的相干结构。图31（a）显示了特征频率St=0.1时的第一SPOD模式。在JL的右分支附近，我们观察到一个从SL末端到尾流区域的波包。图31（a）–31（d）的比较表明，该波包的空间波长和跨度随着St的增长而下降。图31（c）和31（d）显示另一个波包平行于第一结构并靠近中心轴。结合图24中Q准则的可视化，我们可以推断出外波包和内波包都与闭合区域附近的涡旋脱落动力学有关。然而，外部波包与SL末端较大的涡旋结构脱落有关，主要包括细长和倾斜的环形涡旋，而内波包与从闭合区中心部分脱落的较小且各向同性的涡旋结构有关（r/dc<0.5）。我们将这两个相干结构分别称为JLout和JLin，如图31（c）所示。请注意，在St=0.1时，JLout比JLin占主导地位，而在较高的频率范围内，JLout和JLin都有明显的体现。图31（c）显示了与SL平行分布的波包（表示为SLKH）。由于剪切效应在该区域占主导地位，我们可以推断这种结构是由空气-水界面上的剪切引发的开尔文-亥姆霍兹（KH）不稳定性引起的。在之前的实验结果中，已经观察到在通气腔中该区域附近的类似KH诱导波结构（参见例如Wu等，2019）。此外，在腔内，我们观察到不规则分布在JL左分支附近的波包（表示为JLr），如图31（c）所示。这些结构与折返涡旋有关，旋涡是由闭合区的湍流混合产生的，并且由于强大的压力梯度而向上游对流。总之，图31中的结果和SPOD分析得出以下结论。首先，观察到与不同机制相关的相干结构，包括KH不稳定性诱导的SLKH、与大尺度涡旋脱落相关的SLout、与小尺度涡旋脱落相关的SLin以及与更发达的湍流上游对流相关的SLr。其次，在特征频率St=0.1时，JLout在流动中起主导作用。最后，SLKH结果与先前的实验结果一致。本研究对雷诺数为O（103）的通气空化流进行了理想化的DNS，忽略了重力效应和界面张力效应。先前的实验研究，例如Mäkiharju等（2013）的研究，已经表明雷诺数可以影响通气空泡流的动力学。为了理解当前DNS方法中采用的低雷诺数的含义，本节对模型相对于实际流量的保真度进行了批判性评估。解决本研究的局限性对于理解模拟结果在多大程度上代表实际流动物理场至关重要。在他们的综合实验工作中，Mäkiharju等（2013）研究了通过空气喷射通气的壁式边界部分空泡闭合处的空气夹带，研究了从106到107的雷诺数范围内几何相似的流动。他们以两个数量级（103至105）的韦伯数变化进行了Froude尺度实验，观察到腔界面模式，闭合拓扑结构和气体夹带率的显著变化。他们发现，较低的Re情况表现出玻璃腔界面，而较高的Re情况则显示出与界面附近液体湍流的更强相互作用。当前DNS研究的局限性主要来自较低的雷诺数，它可能无法在实际强度下捕获液相和气相之间的复杂相互作用。此外，忽略重力和表面张力效应可能会影响结果对真实世界尺度和几何形状的适用性。因此，有必要在更高的雷诺数下进行进一步研究，并结合其他相关效应，如弗劳德数和韦伯数，以便更全面地了解通气空化流。尽管当前研究中的雷诺数低于实际值，但DNS仍然可以为通气空化流的基本动力学及其与湍流的相互作用提供有价值的见解。例如，一些研究已经检查了雷诺数在破碎波流中的影响，重点关注粘性耗散、穿透深度以及从平面流到三维流的过渡等关键特性。Deane，Stokes和Callaghan（2016）讨论了湍流饱和度，其中羽流横截面积的增加比羽流内流体剪切应力的增加更相关。Giorgio，Pirozzoli和Iafrati（2022）证实了耗散能量分数与雷诺数或维数无关的发现。Mostert, Popinet & Deike（2022）观察到，雷诺数不会影响Re达到阈值后气泡尺寸分布的形状，因为此时平均湍流耗散率对Re不敏感。本研究讨论的通气空泡流动和破浪流都是空气-水湍流混合流系统。尽管它们的物理特性不同，但破浪研究的结果可以为通气空泡研究提供有价值的见解。根据Deane et al. （2016）、Giorgio et al. （2022） 和 Mostert et al. （2022）的研究结果，本研究获得的空气的平均体积分数、平均速度和泡内压力分布可能有助于理解这些流动过程的潜在动力学。此外，本DNS可以通过涡流动力学、雷诺应力、TKE产生和传输、泄气和涡脱落过程以及气体交换、涡流与雷诺应力之间的关系等详细分析，为流动动力学和湍流提供有价值的信息。总之，本研究中的低雷诺数模拟具有揭示通气空化流基本动力学的潜力。认识到所选雷诺数的局限性以及对表面张力和重力影响的忽视，未来对通气空化流的研究必须考虑这些影响。通过考虑当前研究的优势和局限性，未来的研究可以更好地理解这种流动的复杂性质，并继续推进该领域。在这项工作中，我们进行了DNS研究，并相应地评估了与湍流特性相关的离散化，以研究通气空化流。利用CLSVOF方法捕获了气相和水相之间的界面，该方法可以准确描述界面几何形状，显示出良好的质量守恒性能，并稳健地捕获两相流动力学特征。将基本非振荡方案应用于对流项，以提高数值精度。空化器主体使用浸没边界法表示。DNS结果提供了湍流通气空化流的详细描述。我们的分析侧重于湍流和通气空化之间相互作用的基本动力学。基于DNS数据，采用时间和空间平均计算了流体相的平均体积分数和平均速度。由于空气和水的剧烈混合，在空泡的大部分封闭区域，空气的平均体积分数约为0.5。在空泡的前部，流向流动相较其他径向流动更占主导地位。平均流展示了几种特征流动结构，包括剪切层（SL）结构，其中空气被外部水流剪切并向前移动，以及再循环区（RA），其中空气向后流向空化器并减速。此外，在空泡的闭合区域，SL后面的水射流，即射流层（JL），是另一种特征结构，其中射流的一个分支向后移动并重新进入气泡，而另一个分支向前移动进入尾流。通过观察沿中心轴的平均空化数来研究空泡内和闭合区附近的压力分布，该空化数保持约σ=0.15的恒定值，然后在闭合区附近由于JL结构内部的停滞区而急剧减小。空化数和空泡剖面的变化与现有的半经验理论一致。分析了平均旋涡及其成分，阐明了空化流中的涡旋动力学。结果表明：由于空气和水之间的剪切效应，主要由横向涡度分量贡献了该涡度分量，该剪切效应集中在空泡前部的SL和整个封闭区域。前部的流向分量可以忽略不计，但在闭合区域增加。通过研究内旋收支方程，研究了空泡内的涡旋动力学。我们强调评估控制涡度产生的项，即拉伸效应项、气压斜效应项和粘性效应项。在空泡的前部，涡旋动力学受所有三种效应的控制，而在闭合区域，涡旋动力学由气压斜效应和粘性效应主导。此外，对雷诺应力的研究表明，在大部分空泡中，峰值出现在空气-水界面附近的流向分量占主导地位。总体而言，雷诺应力在闭合区域附近变得更强，其中对角线分量具有可比性，表明该区域的流动变得各向同性。此外，我们发现TKE在空泡的闭合区域附近很大，特别是在JL结构中。此外，在SL和JL结构中，流向湍流波动占主导地位。此外，SL和JL结构在TKE的产生和传输中起着不同的作用。在SL结构中，由于剪切效应，TKE在生产传输项上增加，而在JL结构中，由于流向拉伸效应，TKE降低。在SL和JL中，压力传递项往往会增强TKE。此外，平均流量对流输运和压力输运项有利于将容器电容器从SL输送到JL。本文还研究了封闭区域的泄气和涡脱落过程。涡流结构的产生与气泡的泄漏有关。此外，这些过程发生在停滞的高压区域，空气和水相剧烈和湍流混合。此外，还计算了闭合区附近水流速和VOF函数的频谱。闭合区的涡旋脱落呈现半周期性模式，能谱呈现内在特征频率。闭合区控制体中空气体积分数的自相关函数表明，泄气和涡脱落过程具有相同的特征频率。为了进一步研究闭合区附近的结构，并将这些结构与涡脱落和泄气过程联系起来，我们进行了SPOD分析，以识别以不同频率振荡的相干结构。能量的累积分数结果表明，第一个SPOD模式比其余模式对能量的贡献更大，特别是在低频范围内，表明低秩行为。前导模式的相对贡献随着下行距离的增加而增加。此外，不同频率径向速度的第一SPOD模态的等值线显示，各种相干结构在闭合区附近具有不同的机制，包括KH不稳定性诱导的SLKH、大尺度涡旋脱落引起的SLout、小尺度涡旋脱落的SLin、上游湍流对流的SLr等。在特征频率St=0.1附近，JLout结构在动力学中起主导作用。我们注意到，本研究使用不可压缩流动模型，没有考虑可压缩性效应。Ganesh，Mäkiharju和Ceccio（2016）的一项研究表明，气泡冲击的形成和传播可以在空泡的动力学中发挥重要作用。他们的研究表明，体积分数约为0.1的空化流具有较低的声速和相对较高的马赫数。由此产生的气泡激波的形成和传播对空泡的动力学有相当大的影响。此外，气泡冲击会改变闭合模式。他们的研究指出了两种不同的机制，这些机制与不稳定空化流中的空泡脱落有关：RJ诱导的脱落和气泡冲击传播，后者显示出强烈的周期性云脱落行为。他们的结果表明，需要通过结合气泡冲击效应来改进空泡和闭合的动力学模型。此外，Ganesh等（2016）发现，与一维气泡激波传播相关的分析模型和基本跳跃关系可以捕获基本的空泡动力学。在空化流中，存在大范围的湍流运动、空气和水之间的剧烈混合以及不同特征结构之间的相互作用。高保真数值模拟通过提供非侵入式方法来详细描述流场，从而补充了实验研究。这项工作是具有合理实践的DNS。在今后的工作中，为了进行更真实的数值模拟，并改进这些模拟在研究通气空化问题中的作用，应考虑以下主题。首先，在未来的应用中需要考虑高于本研究中使用的雷诺数。其次，我们强调所研究的问题是主动通气，而不是自然空化。未来的模拟应考虑自然空化的影响，以研究其对通气空化流动力学的影响。此外，应开发适当的两相子网格尺度湍流模型，用于空化流的大涡模拟。此外，还需要根据SL和JL结构的特性建立降阶模型，用于超空泡优化等应用。此外，气体交换、闭合涡流和雷诺应力之间的关系需要通过使用更广泛的参数进行更多的模拟来进一步研究。此外，未来的模拟应考虑可压缩性，以研究气泡冲击对空泡的影响（Ganesh 等 2016 ）。本研究中采用的计算方法已由我们小组在先前对其他多流体问题的研究中得到广泛验证（例如 Hu等，2012; Yang 等 2017，2018; Gao 等 2021a,b）。为了进一步验证数值方案在解决当前通气空化问题方面的适用性，我们通过将结果与文献结果进行比较，进行了更多的测试，并进行了网格收敛测试。请注意，尽管我们的模拟通过求解Navier-Stokes方程来考虑粘性效应，但空泡的前半部分几何形状主要由无粘性流动动力学主导，这意味着我们的CFD结果可以与势流理论之间进行以下比较来验证。
图32：给定阻塞率B的预测最小空化数σmin。红色方块、绿色三角形和蓝色圆圈分别对应于使用粗、中和细化网格的仿真结果。黑色虚线对应于Brennen（1969）理论中的公式。
Brennen（1969）基于势流理论研究了阻塞比对最小空化数σmin的影响，其中σmin在大Cq值下接近。在数值测试中，我们使用了三组不同分辨率的网格，即Δxmin/dc=0.08、Δxmin/dc=0.04和Δxmin/dc=0.02。图32显示，当网格细化时，根据我们的模拟计算的σmin值接近Brennen（1969）的预测。结果表明，Δxmin/dc=0.02的网格分辨率足以准确捕获空泡形态，在B<18%时。我们还使用与Cao等（2017）实验相同的参数值进行了模拟，即B=9%，Cq=0.14和Fr=18。图33比较了空泡轮廓。我们的结果与空泡前部的实验剖面一致。
图33：Cao等（2017）实验（红色圆圈）和本模拟（黑色实线）中的空泡轮廓比较。参数设置为B=9%、Cq=0.14和Fr=18
为了验证本研究所用三个网格的收敛性，比较了粗网格M1、中间网格M2和细网格M3上的平均流量和模拟湍流统计。网格M1、M2和M3的详细信息在第2.3节中提供。图34显示，M2和M3仿真之间的等值线φ(平均)=0.5的差异比M1和M2仿真之间的差异小得多，这表明在细化网格时有收敛的趋势。
图34：模拟中φ(平均值)=0.5的等值线与网格M1、M2和M3的比较
图35：多相湍流统计的网格收敛结果：（a）φ(平均),（b）u(平均),（c）u’u’和（d）ϵ
图35绘制了沿流向轴的平均流量和湍流统计数据，包括平均VOF函数φ(平均)、平均流速u(平均)、流湍流波动u’u’和湍流能量耗散率ϵ。一阶统计量（φ(平均)和u(平均)）和二阶统计量（u’u’和ϵ）均通过使用网格网格M1、M2和M3进行仿真定性捕获。使用网格M2和M3的结果之间的差异远小于使用网格M1和M2的结果之间的差异，这再次显示了在细化网格网时收敛的趋势。网格M2和M3之间的结果差异可以忽略不计，表明数值收敛。
为了确认本研究中采用的网格分辨率足以捕捉通气腔封闭区域的湍流动力学，我们计算了湍流涡旋的最小尺度，即Kolmogorov长度尺度，定义为（Kolmogorov 1941）根据Pope（2000）给出的标准，如果局部像元大小Δ与η的比值满足Δ/η<2.1，则可以认为模拟是具有Kolmogorov η量表的DNS。图36根据使用网格网格M1、M2和M3的模拟，绘制了腔体闭合区域的Δ/η。结果表明：M3模拟中的Δ/η比值满足Δ/η<2.1，而M1和M2模拟中的Δ/η比值在闭合区中为Δ/η≥2.1。本文利用在精细网格M3上获得的数据进行分析和讨论。
图36：（a）网格M1、（b）网格M2和（c）网格M3上局部像元大小与Kolmogorov标度Δ/η的比率等值线。黑色虚线对应于φ(平均)=0.5的等值线
To examine further the Reynolds shear stress τ R xr shown in figure 14(d), we perform a quadrant analysis to illustrate the contributions by four different quadrants of τ R xr, which are classified according to the signs of the velocity fluctuations: Q1(+√ρ ux , +√ρ ur ), Q2(−√ρ ux , +√ρ ur ), Q3(−√ρ u x , −√ρ ur ), and Q4(+√ρ u x , −√ρ ur ) (Wallace 2016). Figures 37(a) and 37(b) present the joint probability distribution function P( √ρ ux , √ρ u r ) and the covariance integrand ρuxur P( √ρ ux , √ρ u r ), respectively, which are obtained at the location where τ R xr obtains its peak value. Figure 37(a) shows that the major axis of the elliptically shaped P( √ρux , √ρu r ) is inclined towards Q4 and Q2. Additionally, the covariance integrand shown in figure 37(b) indicates that the Q2 and Q4 quadrants contribute more to the Reynolds shear stress τ R xr than the Q1 and Q3 quadrants.为了进一步检查图14（d）所示的雷诺剪切应力τRxr，我们进行了象限分析，以说明τRxr的四个不同象限的贡献，这些象限根据速度波动的迹象进行分类：Q1(+√ρ ux”,+√ρ ur”)，Q2(−√ρ ux”,+√ρ ur”)，Q3(−√ρ ux”,−√ρ ur”)和Q4(+√ρ ux”,−√ρ ur”)（Wallace 2016）。图37（a）和37（b）分别表示了在τRxr获得峰值的位置获得的联合概率分布函数P(√ρ ux”,√ρ ur”)和协方差积分ρux”ur”P(√ρ ux”,√ρ ur”)。图37（a）显示椭圆形P(√ρ ux”，√ρ ur”)的长轴向Q4和Q2倾斜。此外，图37（b）所示的协方差积分数表明，Q2和Q4象限对雷诺剪切应力τRxr的贡献大于Q1和Q3象限。
图37：（a）联合概率分布函数P(√ρ ux”,√ρ ur”)和（b）协方差积分ρux”ur”P(√ρ ux”,√ρ ur”)在τRxr达到其峰值的位置的等值线。（a）和（b）中的等值线区间分别为[0.01, 0.1]和[0.02, 0.1]。
本文翻译转载自《JFM》杂志“Simulation-based study of low-Reynolds-numberfow around a ventilated cavity”。
作者：Han Liu，Zuoli Xiao，Lian Shen
