船海结构物水气泡混合流机理及数值模拟方法研究进展

张晓嵩，和康健，万德成*,3

1 上海交通大学 船海计算水动力学研究中心，上海 200240

2 上海交通大学 船舶海洋与建筑工程学院，上海 200240

3 浙江大学 海洋学院，浙江 舟山 316021

0 引 言

在海洋环境中，船舶与海洋工程结构物在与周围流体相互作用时，会出现明显的白色泡沫流动，即为船海结构物水气泡混合流。在实际作业过程中，船海结构物水气泡混合流广泛存在且形式多样。对于运动中的船舶，由于船体与流体之间的强烈撞击，船艏部会出现兴波破碎，产生气体卷吸和液体飞溅，从而形成水气泡混合流[1]；在船艉部，由于船艉后方的流动分离，以及螺旋桨空泡流动与湍流的相互作用，会引起剧烈的水气泡多相流的相互渗透与混合，从而形成大范围扩散分布的白色气泡尾迹[2]。船舶周围产生的水气泡混合流主要由船舶的航行运动引起，而海洋平台周围出现的水气泡混合流则是由平台与波浪间剧烈的相互作用引起。例如，波浪砰击平台立柱或平台下甲板时，会导致水气相互掺混与液滴飞溅雾化，形成水气泡混合流[3]；在极端海况作用下，波浪冲击平台主体会产生甲板上浪和剧烈的冲刷流动，从而引起剧烈的水气掺混，形成水气泡混合流[4]。以上现象均为船海结构物水气泡混合流的表现形式。作为船舶与海洋结构物等众多自由液面工程实际问题的重要组成部分，这种水气泡混合流是普遍存在、易于观察的，但也难以系统描述和分析。

从流动组成上看，水气泡混合流现象表现为：在水面以下存在大量的气体卷吸进而演化成多尺度的气泡，水面以上则存在大范围的液体飞溅进而演化成雾化液滴，而气泡和液滴同时在大变形自由面附近发生两相剧烈的掺混，形成水、气、泡、雾耦合的流动状态。

水气泡混合流的生成演化包含着剧烈的多相流动和多物理场的力学特性，对船舶结构物的水动力学性能、极端载荷、流激噪声、尾迹形态等有着重要影响，重点体现在以下方面：

1） 水气泡混合流是船海结构物直接与周围流体剧烈地相互作用而生成，会引起船海结构物能量损失或是载荷变化，直接导致船舶快速性能和操纵性能的降低，或海洋平台结构物耐波性能和稳定性能的改变[5]。

2） 水气泡混合流中的分散多相流动与湍流相互作用，会引发流动剧烈脉动，改变船海结构物周围流激噪声的辐射特性，并会对声学探测设备的信号传输产生影响[6]。

3） 水气泡混合流会导致波浪砰击结构物的载荷特性与摇荡运动特性发生变化，飞溅高度增加，载荷扩散范围增大，并加剧甲板上浪现象[7]。

4） 舰船尾迹中的水气泡混合流容易被探测和追踪，致使舰船隐身性能大大降低[8]。

由此可见，船海结构物水气泡混合流与船海工程装备的绿色节能、减振降噪、结构安全、尾迹隐身等密切相关，是船舶与海洋工程领域需要解决的重要问题。

当前，国内外针对船海结构物水气泡混合流的研究热点和难点问题主要包括2 个方面：

1） 水气泡混合流的生成演化及其与船海结构物相互作用的机理问题。水气泡混合流本身的物理机制高度复杂，涉及多相流、多尺度、剧烈湍流等问题的耦合，并且受多种影响因素的控制。考虑到其与船海结构物的相互作用，还需要探明水气泡混合流与多系统运动边界、绕流和分离流的耦合，这些问题将研究的复杂性提升到了更高水平。当前，机理探索方面的研究可以分为单纯水气泡混合流生成演化机理和船海结构物周围水气泡混合流机理两方面。对于前者，需重点从微观层面分析因多相流动界面撞击作用而产生的气泡和液滴结构，并从物理特性的角度研究水气泡混合流的多因素影响规律；对于后者，则重点从宏观层面分析船海结构物不同位置处水气泡混合流的生成条件与演化分布，从流动、能量和载荷等方面研究水气泡混合流与结构物性能的关系。机理问题的研究有助于深化人们对船海结构物水气泡混合流问题的认识和理解，并为进一步对其进行控制和优化奠定理论基础。

2） 水气泡混合流的高可靠性数值模拟方法问题。水气泡混合流的一个重要特点是其广泛存在于实尺度结构物周围，而在模型尺度试验中则难以复现。在相似规律不明确的条件下，难以通过试验手段对实际情况进行预报，制约了在船海结构物水气泡混合流方面的研究突破。因此，发展高可靠性数值模拟方法就显得尤为重要。当前，在数值方法方面的研究可以分为精细化模拟方法和模式化模拟方法。对于前者，重点是建立可高精度捕捉微尺度界面变形和气泡、液滴结构生成演化的数值算法，解决界面捕捉、数值耗散、守恒格式等关键问题，用以进行船海结构物水气泡混合流的微观物理机制分析；对于后者，重点是建立能模拟实际大规模气泡液滴群的精确物理模型，并与界面流动模拟方法相结合，形成跨尺度流动求解方法，用以对实际船海结构物周围水气泡混合流的宏观特征进行预报。数值模拟方法的研究可以为流动机理分析提供有效手段，也能为开展实际结构物周围水气泡混合流动的预报和构型优化设计提供分析工具。

因此，本文将分别对机理探索研究和数值方法研究的进展予以介绍。

1 船海水气泡混合流机理

1.1 水气泡混合流生成演化机理

水气泡混合流中包含着多尺度的流动结构[6]。在微观上，有毫米及亚毫米尺度的微气泡和微液滴；在宏观上，有大气泡、液体射流和大范围的水气自由面。微观尺度的流动结构在表面张力的作用下可以保持一定的形状，而宏观尺度流动结构的特点是水气界面存在连续不断的变形。在生成演化机理探索上，目前国内外的研究可以分为4 个方面：一是以宏观流动为基础，研究局部区域内相界面的撞击作用和破碎过程，并分析小尺度流动结构的来源；二是以微观流动为基础，研究气泡、液滴结构的撞击作用和演化形式；三是研究黏性和表面张力对水气泡混合流的影响；四是研究水气泡混合流与噪声的相互作用。下面，将从这4 个方面分别介绍相关的研究进展。

1.1.1 局部区域多相界面撞击作用

船海结构物周围局部水体的水气界面撞击作用，是船海结构物水气泡混合流生成的一个重要机理问题。在撞击过程中，由界面破碎引起的非线性作用会产生气泡掺混和液滴飞溅，过去，通常是利用波浪翻卷破碎试验来模拟这种界面撞击现象。在试验中，宏观的波浪破碎流动过程主要采用高速摄像系统进行记录[9-11]。在早期研究中，通过现象观察以及气泡尺寸测量，发现水下气泡群呈现明显的多尺度特性[12]，如图1 所示。气泡直径以Hinze 尺度为界，具有不同的生成机制，较大气泡的生成机制属湍流诱导破碎，其决定了水气泡混合流中大尺度气泡的体积分布；较小气泡的生成机制主要是射流和液滴与自由面的撞击作用。Kimmoun 和Branger[13]采用粒子图像测速（particle image velocimetry，PIV）技术，通过测量水下速度场和漩涡流动，研究阐述了湍流是大气泡破碎的主要原因，以及破碎后不同尺度的气泡所导致的分布区别。直径较小的气泡受到的浮力较小，这些气泡会在海洋中聚集与漂浮，形成气泡层[9]。在自由面上方，波浪翻卷撞击会在水面产生高速射流，并在空气中分裂成大量的小水滴，然后在剧烈的流动条件下，水滴会在空气中雾化，形成雾化水滴云。雾化水滴云和水下气泡云是船海结构物水气泡混合流的重要流动结构。

图1 波浪破碎诱导的水气泡混合流中的多尺度特性[12]Fig. 1 Multi-scale characteristics of water-air-bubble mixed flow induced by wave breaking[12]

由波浪翻卷撞击破碎试验，可获得丰富的水气泡混合流相关的流动现象，但试验中最大的困难在于对微观流动结构无干扰的定量统计和测量。为更加准确地探究波浪破碎的演化及破碎过程中气泡、液滴运动的动力学特性，数值模拟方法被应用到波浪破碎问题研究中。早期，因计算条件所限，部分学者在湍流模拟上采用雷诺平均（RANS）或大涡模拟（LES）等模型化方法进行了模拟计算[14-16]，并结合了一些用于描述离散气泡群效应的模型，重点是分析试验中难以测量的流动速度场，分析波浪翻卷到波前产生飞溅的过程，并提取飞溅液滴以及气泡掺混的情况，但该模型化方法的缺点是流场中的脉动并不能被准确求解，导致气泡多尺度结构特征描述存在一定的误差。

近年来，随着计算技术的快速发展，许多学者[5,17-18]开展了波浪撞击破碎过程的三维直接数值模拟（direct numerical simulation，DNS），典型的模拟结果如图2 所示。从中可以看到，能精细捕捉到撞击过程中的复杂界面破碎。借助DNS，还可以发现波浪破碎早期的涡丝现象，这也进一步解释了波浪破碎过程中的卷气机制以及相应的能量损耗机理。另一方面，DNS 的精细效果也推动研究者可以从宏观−微观的耦合层面来理解微气泡等流动结构的生成。Chan 等[19]的研究表明，在一个27 cm 范围内的水气泡混合流中，每秒可发生上百次的微气泡撞击，可见，在实际中每立方米空间内每秒会产生上百万个微气泡。这种多尺度流动现象的探索对理解实际船海水气泡混合流生成演化机理具有启示作用。

图2 波浪破碎中水气泡混合流生成与涡丝的数值模拟结果[17]Fig. 2 Numerical simulation results of water-air-bubble mixed flow generation and vortex filament in wave breaking[17]

基于多尺度的精细化数值模拟，Wang 等[18,20]总结出了波浪破碎问题中的4 种水下气泡产生机理，如图3 所示，包括射流碰撞卷气、碰撞形成二次飞溅卷气、反向飞溅与初始波峰之间夹带气体（反向飞溅卷气）和湍流夹带卷气。在不同的机理作用下会形成不同尺度的气泡，而这些气泡又会进一步呈现出不同的动力学行为。同时，水面以上雾化液滴的形成也被描述为液体表面拉伸、飞溅凸起和液体脱落。Wang 等[18,20]研究发现，微观流动之间的相互作用发生在自由面附近，即气泡在自由面附近破碎会生成新的液滴，而液滴下落砰击自由面也会产生新的气泡，这是除大尺度流动破碎演化之外微尺度流动结构新的来源。

图3 破碎波水气泡混合流中的4 种卷气机制[18]Fig. 3 Four air entrainment mechanisms for water-air-bubble mixed flow of wave breaking[18]

1.1.2 微观流动界面演化机理

船海结构物水气泡混合流多相流动界面的相互撞击过程实际包含了大量微观流动结构的界面演化，如气泡、液滴等。目前，一些研究人员主要从少量气泡的相互作用，以及气泡和液滴与自由液面相互作用的角度，对气泡液滴生成微观机理进行了研究分析，例如单个大气泡的破碎过程试验[21]、极少量气泡间非线性相互作用试验[10-11]、气泡上升与自由面相互作用试验[22]、液滴下落撞击自由液面过程试验[23]等。这些试验研究利用高速摄相机记录了水气界面复杂演化的全过程，研究分析了气泡破碎过程的发展规律以及末期射流的产生原因、气泡间撞击及聚并过程的非线性作用机理（图4），为气泡云及水气泡混合流生成机制的研究奠定了基础。

图4 微观界面撞击形态关于Oh 数和We 数的相图[23]Fig. 4 Regime maps for micro-scale interface morphology based on Oh and We[23]

在以上试验中，测量设备和仪器会对气泡、液滴的形成产生干扰。为克服该问题，近年来，研究人员开始采用计算流体力学（CFD）数值模拟方法开展研究。在气泡的微观界面撞击作用机理方面，模拟研究了气泡与气泡之间的靠近追逐[24]、碰撞反弹[25]、接触聚并[26]等过程。在靠近追逐问题上，共轴线气泡上升是机理研究的经典问题，当前的研究已经涉及到了各种雷诺数范围[27]，分析了上部气泡的屏蔽作用，并在雷诺数Re≈ 100的情况下发现了稳定“平衡距离”现象[28]。然后以此为基础，又进一步将研究拓展到了并排上升和多角度上升方面，其中尾流效应仍是研究的重点[29]。在接触聚并问题上，目前有已得到广泛认可的气泡/液滴聚并过程描述[30]。对于气泡来说，首先两气泡彼此靠近，在二者之间形成一层非常薄的液体层，然后随着压缩作用液体的排出，在达到临界厚度时破碎，从而产生局部的微观气桥现象，最后两气泡中的气体相互连通，即完成聚并。对于液滴，也有相似的气体层和液桥现象[31]，主导着聚并过程的演化。聚并转化的条件受表面张力，也即韦伯数We的控制。

在当前基于网格类方法的数值模拟中，大多是采用流体体积（volume of fluid，VOF）法或是水平集（Level-set）法进行水气界面的捕捉，但这两种方法均存在一定的不足。VOF 法采用连续流体模型计算界面表面张力存在较大误差，而Level-set 法在保证质量守恒上存在问题。因此，有学者耦合这两种方法的优势，采用CLSVOF 方法进行了比较复杂的气泡/液滴微观多相流动界面撞击流动模拟，包括单液滴撞击倾斜液膜飞溅过程模拟[32]、三气泡上升过程中的相互作用模拟[33]等。借助于该方法，可以更好地描述气泡/液滴表面的变形和运动特性，从而能为进一步开展表面张力对流动特性的影响研究奠定基础。另一方面，也有研究采用粒子法对高度分散的微尺度气泡液滴生成演化过程进行了模拟[34-38]，并对翻卷破碎、射流等剧烈的自由面演化以及由其引起的流场脉动问题进行了分析。

1.1.3 黏性和表面张力影响机理

水气泡混合流的生成演化受流体黏性、湍流流动以及界面之间表面张力等因素的影响，在有关影响因素的研究中，由于其参数调整的便捷性，理论研究和数值模拟成为主要的研究手段。同时，也有一部分试验工作是通过采用大部分物理参数接近而某一物理参数不同的流体介质来对影响机理进行研究。

在有关黏性影响效果的研究中，Zhang 等[39-40]以不考虑黏性的势流理论为基础，通过引入黏性修正压力和正应力的黏性分量对非零黏性剪应力进行补偿，分析了有、无黏性导致的区别。该研究分析了单个气泡在浮力作用下的变形、破碎以及射流过程，发现在考虑黏性的情况下，气泡射流速度是随流体黏性的增加而衰减的，气泡射流能量降低，流场损失的动能和势能逐渐转化为黏性耗散能。同时，还研究了黏性对单个气泡和双气泡脉动特性的影响规律，发现流体黏性会抑制气泡脉动。以上研究的特点是计算效率比较高，可以开展系列研究并得到宏观的规律。王昭太等[41]则采用基于黏流CFD 的模拟方法，探讨了流体黏性对气泡上升运动特性的影响，研究认为流体黏性主要影响气泡内的速度分布，阻力系数随黏度比的增大而减小。除了针对单个气泡/液滴动力学行为的研究，有学者还通过数值模拟分析了流体黏性对水跃（hydraulic jump）等复杂界面流动中水气泡混合流气泡掺混特性的影响[42-43]，结果表明，随着黏性作用的减小，掺混区内最大空隙率和气泡数密度呈单调增加的趋势。黏性影响和湍流效应密不可分，Thoroddsen 等[44]利用PIV 技术，开展了湍流流动中气泡运动变形的力学行为研究，并得到了指定区域的剪切速率、湍流动能和动能耗散等信息，发现气泡在湍流场中的破碎演化过程是由外部因素（如湍流强度、剪切速率）和内部因素（连续相的流体动力学因素和相界面相互作用）共同作用的结果。Zhang 等[45-46]通过数值模拟分析了湍流边界层与气泡群的相互作用，并给出了湍流作用下的气泡迁移运动轨迹与尺寸分布特性。

对于表面张力影响研究，在单独气泡/液滴动力学行为方面，艾旭鹏和倪宝玉[40]研究了表面张力（We数）对单个球状脉动气泡的影响，结果显示表面张力会缩短气泡的振荡周期但幅值不变， 同时还会使势能增加。Ma 等[47]研究了表面张力对从边界注入的单个气泡生成过程的影响，结果显示表面张力的增加会导致气泡生成平均周期增加、气泡体积增加、气泡生成频率降低。在气泡群动力学行为方面，Gourich 等[48-49]以气泡羽流问题为基础，分析了流体表面张力对气含率、气泡形态和尺寸分布的影响，认为表面张力的变化会明显影响气泡群中的聚合现象，气泡聚并速率将随表面张力的降低而降低，水、气两相体系的传质速率也会明显降低，并且流场中的气泡体积分布更小。Chaumat 等[50]进一步研究了表面张力对气体破碎的影响，发现表面张力减小后，会使水气泡混合流在低流动速度下更容易发生失稳，进入液体流场中的气体破碎得更快，这同样也是气泡尺寸分布减小的重要原因。对以上文献中介绍的影响因素、研究对象和影响效果予以了总结，结果如表1 所示。

表1 黏性和表面张力对水气泡混合流的影响效果Table 1 Effect of viscosity and surface tension on water-airbubble mixed flow

1.1.4 对噪声的作用机理

水气泡混合流动内部包含数目庞大的气泡和液滴微尺度流动结构，致使整个流场存在高度的脉动。这一特点使得水气泡混合流会对水面船舶噪声的产生与辐射产生重要影响。

首先，水气泡混合流中气泡和液滴的动力学行为会对流场脉动造成影响，产生新的声源。Song 等[51-52]利用水听器对气泡形态变化过程中周围声压的变化进行检测，并在试验中分析了气泡大小、直径、溃灭方向等因素对噪声辐射的影响。同时，也有学者认为液滴产生的噪声本质上也来源于气泡[53]。当液滴撞击自由面时，引起的掺混气泡受到局部压力的作用会产生振荡，进而辐射噪声。在声理论中，考虑将气泡构造为声源的方法可以分为线性理论和非线性理论。当前，主要的方法是线性理论，即只考虑球状收缩膨胀的单极子源[54]。但近年来的实验研究证明，对于复杂的流动，非线性声源对水下噪声的影响同样不可忽视。

其次，气泡的存在会对噪声的传播产生影响，这主要体现在声速变化和声功率损失2 个方面[55-56]。当噪声在流体中传播穿过气泡时，气泡会影响声速并耗散能量[57]。早期，Commander 等[58]针对线性声波穿过水气泡的问题进行研究，提出了一套预测方法。 在该线性波理论的基础上，Louisnard[59]开始关注非线性成分对声波传递的影响，并研究发现，在气泡流中，非线性波相比线性声波衰减得更明显。随后，Fuster 等[60]在此基础上又进一步考虑了质量转换、气泡融合等动力学行为，发现可以更精确地预测水气混合流环境下的噪声。

1.2 结构物周围水气泡混合流机理

船舶和海洋平台作为典型的船海工程结构物，在实际作业过程中，其周围广泛存在着水气泡混合流现象，不过机理和作用效果有所不同。船舶周围的水气泡混合流主要由其航行运动引起，并具有明显的位置区域差异性，即在船艏部、舯部和艉部不同位置处生成的水气泡混合流现象与机理均不同。而海洋平台周围的水气泡混合流则主要源于恶劣海况下的波浪砰击和甲板上浪。下面，将从这几方面来介绍相关的研究进展。

1.2.1 船艏水气泡混合流现象机理

船艏兴波是伴随船舶航行运动而产生的，在实尺度下一般不是规则理想的波动，而是会出现船体艏部与流体的剧烈撞击，出现波浪破碎和水气渗混，进而形成水气泡混合流。因此，船艏部的兴波破碎成为研究的首要问题。经典的试验研究主要针对的是不同弗劳德数Fr下的船艏兴波形态[61-63]。在试验中，Fr的范围一般在0.28～0.45之间，可通过提高航速来诱导产生兴波破碎，然后采用自由液面可视化技术观测船艏兴波破碎产生的过程、艏波破碎流动结构以及艏波在船体表面上的分离点，还可以利用PIV 技术测量船模周围的速度场和涡量场，分析在兴波起始位置处出现的负涡区和在船边界波峰处出现的正涡区。试验研究指出，船艏的几何构型特征会使兴波向两侧翻卷，在高速情况下，兴波的剧烈翻卷撞击是产生船艏部水气泡混合流的主要原因。

前期的实验研究针对多种船艏几何构型和航速，给出了兴波破碎生成的水气泡混合流现象的结果，再运用数值模拟方法，可以进一步定量研究船艏兴波破碎产生的水气泡混合流的物理量。由于需要考虑复杂几何曲面结构物的影响，数值模拟研究的难度要比前文所介绍的纯波浪破碎模拟大得多，因此前期主要的数值模拟研究在湍流模拟方法上是以RANS 为基础[64-65]，Fr数最高可达0.62。通过模拟，可以成功复现兴波在船艏侧面的翻卷和射流现象，但平均化方法对于破碎的细节却难以捕捉。近年来，基于分离涡模拟（DES）的湍流模拟方法成为船艏水气泡混合流模拟的主要手段[66-69]，Fr数最高进一步提高到了1.2。典型的数值模拟结果与实船海试结果的对比如图5 所示，图中A，B，C，D 分别为4 个典型的现象对比区域。从中可以看出，由数值模拟可以得到船艏射流分离、自由面大变形、兴波多次翻卷以及水气掺混现象。这些数值研究是利用轴向涡量、速度和自由面表面拓扑的截面分析来解释兴波破碎的演化的。首先，液体在船艏表面爬高并向外侧翻卷形成一次破碎，在高速情况下，射流在与自由面接触后会发生二次甚至是多次反射，形成多次破碎，并在射流与自由面接触的位置留下明显的痕迹。通过对涡量场的分析，发现射流末端以不连续的速度与自由面重新接触后，在该位置会形成一个强反旋转涡对，从而进一步形成了大横流梯度和新的射流。从现象完整性的角度来看，数值模拟研究可以用于分析讨论船艏水气泡混合流生成前的过程和机理，但还不能对破碎后形成的大范围水气泡混合流进行准确复现。另一方面，当前的试验和数值模拟强调在高Fr数下才能出现破碎现象，但实尺度船舶在以低速航行时仍会出现水气泡混合流，这说明该现象并不能完全满足弗汝德数相似，其深层的影响因素和流动机理仍需进一步研究。

图5 船艏兴波破碎水气泡混合流数值模拟与实船海试对比[67]Fig. 5 Comparison between numerical simulation and sea trial of water-air-bubble mixed flow induced by bow wave breaking[67]

波浪与船艏的砰击作用是船艏生成水气泡混合流的另一个重要原因。在波浪的作用下，船艏抬升和下落时与水体撞击会产生强烈的砰击雾化和水下气泡掺混，从而形成大范围的水气泡混合流。近年来，Bodaghkhani 等[70-72]利用先进的测量设备针对该现象开展了模型试验研究。在试验中，基于规则几何平板与复杂几何船舶模型，对波浪砰击过程进行了研究，并采用气泡图像测速仪（BIV）和数字粒子图像测速法（DPIV）这2 种图像测速技术研究了破碎波在冲击平板的过程中产生的雾化现象。值得注意的是，该研究结合了实船海试中的一些流动现象，通过与模型试验相配合，对砰击雾化现象建立了新的认识。针对船艏砰击形成雾化液滴的问题，已形成了一些有效的理论模型[71]，并已针对雾化液滴尺寸、速度和空间分布进行了建模，可用于流动现象的快速预报。

1.2.2 船舯水气泡混合流现象机理

船舶舯部附近是船舶水气泡混合流演化的重点区域，在船艏形成的水气泡混合流会在船舶舯部附近形成一定的分布。前面已提到，在模型尺度下由船模试验产生水气泡混合流需要很高的航速，且生成的水气泡混合流也很难演化覆盖整个船体区域，因此难以从模型尺度的角度研究该演化问题。Johansen 等[1]开展了很有价值的实船观测研究，在研究中对实尺度下“雅典娜II”号船周围的气泡流场进行了测量，提出了一种考虑实际探头尺寸的气泡测量新方法，并通过实测对流场中的气体体积分数、气泡速度和气泡尺寸分布予以了统计。经研究发现，由于船体周围的加速度，在船舯处演化的气泡平均速度要略大于船速，并在边界层内部逐渐减小；气泡速度随着深度和侧面距离的增加而增加，并且靠近自由面附近气泡的平均尺寸更大；当在不同航速、不同位置逆浪航行时，掺气量均随航速的增大而增大。另外，研究还针对湍流边界层气泡卷吸以及尾部气泡流动等现象进行讨论，得到了气泡速度关于弦长的分布，称为“弦长分布”。

Ma 等[72]基于实船观测，提出了Phenomenological Air Entrainment 模型，用于预报水面船舶周围的气泡，得到了船体周围气相分数分布规律，如图6（图中，εg为气泡源项强度）所示。同时，还讨论了在直航和定常回转运动中，夹带气流的位置与产生气泡流间的定性差异，分析了船舶附体周围的气泡聚集现象，并将仿真结果与海上实验观测结果进行了比较，所得结论可为后续通过数值模拟进行机理探索奠定了基础。船舶舯部附近的水气泡混合流中的气泡掺混量极大，采用界面捕捉类方法模拟计算的成本很高，因此对于船舶周围气泡流的数值模拟研究，主要采用以双流体方法（two fluid model）为基础的模型化方法[73]，以便于分析船舯水气泡混合流的宏观分布和演化过程。

图6 船舯水下掺混气的数值模拟结果[72]Fig. 6 Numerical simulation results of air entrainment at the midship positions[72]

1.2.3 船艉水气泡混合流现象机理

船舶艉部水气泡混合流最重要的宏观表现是船后大范围的白色泡沫尾迹，且该处的水气泡混合流同时还受船舶尾部湍流和螺旋桨扰动的影响。由于该问题对于船舶的多种水动力性能和声光特性具有重要影响，近年来得到了重点关注。对于船尾白色泡沫尾迹的生成机制，目前的研究认为其主要有2 个方面来源：一是航行过程中由船尾附近自由面流场翻卷引起的气体卷吸[74]；二是由螺旋桨高速旋转所产生的空化气体[75]。自由面卷吸的气体将以气腔的形式进入水下，同时还会受到海水压力和船尾分离漩涡流动的影响，导致其表面张力因难以克服流场作用力而发生破碎，最终成为微尺度气泡。空化的气泡则源自水中固有气核在压力低于饱和蒸气压条件下时的汽化，同样表现为微尺度气泡。这两方面的气泡会在船尾形成具有一定长度的气幕带，也即白色的泡沫尾迹[76]。

在实测研究方面，重点关注了大尺度结构和宏观参数。Peltzer 等[77]在早期对高速船尾部的水气泡混合流进行了海上观测，研究提取了船舶尾流中不同位置的多个横向剖面进行流场分析，发现在所有尾迹剖面中均存在明显受到表面张力作用的带状流动结构，该特征在船尾后方数十公里处仍然存在。通过声呐技术，可以对尾迹水气泡混合流进行成像[78-79]，Soloviev 等[78]利用三维声呐技术，以内河航道中的船舶为基础开展了实船试验，分析了因自由面破碎、螺旋桨扰动所产生水气泡混合流现象间的区别，进一步阐述了船体周围气体卷吸的位置和过程，并将实船尾迹的流动结构分析精细化到厘米尺度。

船舶尾部的水气泡混合流现象与船艉型线有较强的相关性。许多试验研究通过高速方尾模型试验研究了尾流的特性。Shen 等[80]利用DPIV 方法对几种不同几何形状船模的尾流水气泡混合流进行了测量研究，重点关注了自由面湍流问题，并在研究中通过数据处理得到了船尾后方自由面附近流动结构和湍流的统计特征，包括不同截面位置的平均速度和脉动速度场。同时，还开展了基于Orr-Sommerfeld 方程的流体稳定性理论分析，发现在实验中观察到的弯曲特征是由平均尾迹的固有不稳定性所引起，导致流动出现了非线性特征。美国海军研究部基于一个方尾的大尺度船模开展了不同航速下的高速拖航试验[81-82]，其利用高速摄影、水下摄影、孔隙率探针、激光雷达等先进测量设备，从含气率、气泡流动、数密度变化、时空分布等角度研究了水气泡混合流的特性，并分析了在不同船速下的区别和影响，典型的试验图像如图7 所示。

图7 高速方尾船模试验中的水气泡混合流[82]Fig. 7 Water-air-bubble mixed flow in the experiment of highspeed square-tailed ship[82]

在数值模拟研究方面，目前已有模型化方法的应用，主要以RANS 或LES 湍流模拟方法与双流体两相流模拟方法的结合为主[83-84]。Castro 等[84]重点讨论了尾流不同位置处气泡体积分布、湍流与气泡流的宏观相互作用和不同航速的影响规律等，结果显示气泡体积与尾流流向位置之间的关系符合指数函数，与船尾相近位置处的气泡体积更小。另外，研究还将尾流中水气泡混合流的生成演化区分为了几个独立的过程，包括由自由面附近由单涡引起的气泡卷吸、复数涡对气泡卷吸的叠加作用、气泡群尺寸的分布及水下扩散过程，完成了物理建模并与试验结果进行了对比验证。Hendrickson 等[85-86]利用隐式大涡模拟（ILES）和浸入边界法（IBM）对方尾尾流水气泡混合流进行了精细化的模拟，讨论了半宽吃水比对尾流演化形态的影响规律，通过气泡识别算法提取了不同情况下水气泡混合流的气泡尺寸、含气率和分布规律，并以水跃问题为参照研究了速度场、两相流分布和湍流特征间的关系。

1.2.4 平台水气泡混合流现象机理

海洋平台因长期作业于固定海域上，容易受到恶劣海况的影响，而波浪砰击和甲板上浪是海洋平台周围水气泡混合流生成的主要原因。目前，相关的机理研究主要针对气泡掺混对结构物所受流体载荷的影响。在大波陡畸形波浪砰击海洋工程结构物的过程中，有一部分气体会被卷入水中，形成多尺度气泡掺混的水气泡混合流，这些气泡会导致脉动增强，从而进一步产生周期性的载荷。图8 所示为波浪砰击直墙结构的载荷测量结果[87]，由图中可以明显看出水气泡混合流的影响，而采用传统理论模型（蓝色曲线）则无法预报该现象。Hu 等[7]在实验中制造出了不同程度的破碎波，用以研究不同状态水气泡混合流砰击海洋平台壁面时产生的载荷特性，包括轻微破碎掺混、带射流翻卷、含大尺度气腔和完全破碎掺混4 种情况。结果表明，含有大尺度气腔的砰击过程在直墙结构上产生的载荷最大；在轻微破碎掺混情况下，波浪爬升的高度最高。同时，因存在气泡掺混的情况，给载荷带来了明显的脉动。该研究说明了不同尺度气泡在水气泡混合流砰击过程中的作用与区别。Sun 等[88]针对掺混气泡波浪砰击固壁表面的问题建立了边界元分析模型，其中重点考虑了翻卷中心的大气腔影响，忽略了边界破碎的小气泡影响，同时分析了气腔体积对压力的影响。目前，大多数相关问题的研究还是以考虑大气腔为主，而采用模型试验和数值计算均很难复现出实际高度破碎的多尺度水气泡混合流现象[70]，因此对于其影响机理还需要进一步开展方法和理论研究。

图8 气泡掺混引起的砰击载荷变化[87]Fig. 8 Effect of air entrainment on impact load[87]

除了由于波浪砰击产生的水气泡混合流外，有研究人员还重点关注了月池结构处产生的水气泡混合流[72]。在当前大型探测船和海洋平台结构物上的大尺寸月池中，自由面的演化过程十分复杂，产生的水气泡混合流会在水下扩散，进而干扰水下设备的使用。Yao 等[89]针对该问题提出了一种抑制气泡群产生和下扫的阻尼装置，并且通过实验分析了水气泡混合流的流动特性和装置的抑制效果。

2 水气泡混合流数值方法

2.1 水气泡混合流精细化模拟方法

船海水气泡混合流现象中包含着复杂的相界面变形和气泡、液滴生成演化物理机制问题，同时，还涉及到多尺度两相流和黏性流、湍流、分离流的耦合作用[18]。为了探明其机理机制，有必要开展精细化数值模拟研究。在网格类方法中，DNS 方法是一种运用三维非稳态的N-S 方程对湍流进行直接数值计算的方法，它不需要引入任何湍流模型，可与高精度的相界面捕捉方法相结合，在研究复杂的两相流动的基本物理过程中起着十分重要的作用。同时，其还能实现对微观流动机制的分析，并对气泡、液滴、湍流涡等难以解析的流动结构进行精准捕捉。DNS 模拟方法对数值格式分辨能力和网格分辨率的要求很高，在算法构建上涉及求解精度、格式耗散、壁面处理、计算效率等多方面的因素，接下来，将从自适应直角网格、相界面求解方法、大密度比两相流求解方法、表面张力模型和浸入边界法5 个核心技术方面介绍相关的研究进展。另一方面，无网格粒子法可以通过自由面粒子来表征飞溅、射流、气泡、液滴等强非线性两相流结构，由于不存在对流项离散耗散等问题，也可以作为水气泡混合流的精细化模拟方法，下面将予以介绍。

2.1.1 自适应加密直角网格

采用传统均匀直角网格开展两相流问题的DNS 数值模拟，在追踪小尺度非定常流动结构时对网格量的需求明显增加，存在计算效率低的问题，为此，发展了自适应网格加密（adaptive mesh refinement，AMR）方法。AMR 通过采用一系列矩形（多层）网格嵌套组合来控制方程在离散求解过程中网格拓扑结构的动态变化，然后根据流场的特点及计算需要对局部区域内网格进行加密，故采用AMR 方法能够在保证在局部网格分辨率较高的情况下提高整体计算效率。当前，已得到广泛应用的AMR 方法主要包含以下3 种：

1） 基于网格点的自适应加密策略，如图9 所示。该技术是以网格点作为树结构的单元，然后针对每个网格点进行划分。Popinet[90]所开发的软件包Gerris 就是典型的基于网格点的AMR 技术，适用于复杂几何外形和界面流动的数值模拟。这款软件包最开始是被应用于具有非规则边界的区域，例如湖泊的浅水方程数值模拟[91]。2009 年，Popinet[92]在Gerris 软件包中植入了VOF 方法和连续表面张力（CSF）模型用来模拟多相界面流动，其利用基于网格点的AMR 网格，得到了一些高精度的结果，如液滴飞溅[93]、波浪破碎问题[11]等。该方法的优点是可针对每个网格单元进行细分，非常灵活，加密效率较高，适用于有限体积方法（FVM）；缺点是必须从不同空间分辨率的相邻网格中取得数据，并行效率低，负载平衡困难。

图9 基于网格点的自适应加密策略Fig. 9 Adaptive grid refinement method based on grid points

2） 基于网格块的自适应加密策略，如图10所示。该技术是以网格块作为树结构的单元，初始计算时，由一个或者多个网格块组成整个计算区域，每个网格块在每个方向上具有相同的固定数目的网格，且每个网格块可以单独进行求解。当1 个网格块被加密时，它会被4 个（二维）或8 个（三维）子网格块所替代，每个子网格块具有和父网格块相同数目的网格且子网格块网格的大小只有父网格块的一半。天体物理的求解程序Flash[94]是一种典型的应用块自适应网格的加密技术。许多研究人员利用基于块结构的AMR 技术求解了多相界面流动问题，包括射流雾化[95]、移动接触线问题[96]等。基于网格块的 AMR 技术的优点是并行效率高、负载平衡容易实现，并且由于每个网格块都是由几何结构完全一致的均匀网格组成，故易于实施高阶差分方法；相对的缺点是其加密操作需要以块为基本单位，计算域整体的加密效率较低，网格量较大。

图10 基于网格块的自适应加密策略Fig. 10 Adaptive grid refinement method based on grid blocks

3） 基于网格片的自适应加密策略[97]，如图11所示。该方法是以网格片作为树结构的单元，既不要求父网格片与子网格片覆盖的范围相同，也不要求每个网格片包含的格点数目相同。该方法综合了前2 种方法的优点，计算区域中的每一个网格都需要单独被检查和操作，加密效率高，但缺点是网格片覆盖的范围变化不定，算法复杂，且网格周围用于插值模板点的数据结构复杂，通信频繁不利于程序并行化，负载平衡较难。

图11 基于网格片的自适应加密策略Fig. 11 Adaptive grid refinement method based on grid patches

2.1.2 相界面求解方法

构造追踪相界面运动的守恒格式，实现高精度的界面重构，是精确模拟水气泡混合流的重要问题之一。采用数值模拟方法实现界面间断关系的先决条件是准确获得相界面的实时位置和形状信息。由于相界面的位置和形状是随时间变化的，因此需要针对界面演化建立模型方程。如果在流动过程中两相流体界面处不发生相变、渗透及相互溶解等现象，可以认为相界面是一个物质面，即界面始终由同一批质点组成。近年来，追踪界面演化的数值模拟方法得到了快速发展，具有代表性的界面数值方法有阵面追踪法（front tracking）[98]、Level-set 法[99-100]、VOF 法[101]和相场法（phase field）[102]。依据相界面的描述方式，界面数值模拟方法可以划分为界面追踪方法（interface tracking）和界面捕捉方法（interface capturing）两类。

界面追踪方法是直接利用物质面假设，采用拉格朗日粒子来描述界面。若要通过当地流场速度推进粒子在流场中运动，可以根据粒子轨迹得到界面随时间的演化情况。阵面追踪法属于典型的界面追踪方法，Unverdi 等[103]将阵面捕捉方法与界面追踪技术耦合得到了一种新的阵面追踪方法，该方法的原理为：一套网格用于求解流体流动；一组自适应元素用来标记运动界面。这样，粒子运动的速度就可以通过欧拉网格的流体速度插值得到。由于所使用的拉格朗日粒子之间的间距一般比欧拉网格的更小，故其界面分辨率较高，可以准确捕捉相界面，但在气泡运动中不能捕捉到气泡的大变形和破碎，以致不能得到广泛的应用。随着粒子法的不断发展，Leung 等[104]使用均匀无网格粒子对界面进行追踪，将粒子与最近的网格进行了关联，结果显示采用该方法可以很好地处理/控制界面拓扑变化。

界面捕捉方法是先定义一个标量场作为相界面的描述函数，然后再采用该标量场的特定等势线来描述界面位置。由此，界面追踪问题即转化为欧拉网格上标量场随时间的演化问题。典型的界面捕捉方法包括Level-set 法、VOF 法和相场法，它们分别选取符号距离函数、流体体积分数和相场函数来作为相界面函数。在标量场时间推进过程中，其局部分布可能会严重偏离理想的分布，从而导致界面法线方向及界面曲率等一些重要界面信息的计算产生误差甚至是错误，因此需要对标量场进行人工修正。不同的界面捕捉方法采取不同的措施来实现：Level-set 法采用重新初始化（reinitialization）方法[105]，VOF 法采用界面重构（interface reconstruction）方法[106]，相场法采用在演化方程中添加包含四阶导数的小量函数方法[107]。VOF 法是一种比较成熟的界面捕捉方法，被广泛应用于流体力学研究与实际工程领域，VOF 分数在网格内的分布示意图如图12（a）所示。基于几何重构的VOF 法具备较好的质量守恒特性，但要实现三维相界面的高阶几何重构比较复杂，目前只有少量针对二维问题的基于最小二乘法的界面重构方法及二次样条函数重构方法[108]。THINC（tangent of hyperbola interface capturing）方法[109]是一种新型的代数VOF 方法，无需界面重构，其特征相函数为双曲正切函数，具体如式（1）所示。由图12（图中，i为网格点序号，C为网格中的体积分数，χ 为相函数）中的对比可以明显看出THINC方法与一般VOF 法在界面表达上的区别，其能有效抑制由界面不连续引起的数值伪振荡，通过结合权重线界面重构方法[110]，能够达到与几何VOF法相当的精度，且算法简单，易于扩展至三维问题模拟中。

图12 网格内相函数分布一维示意图[110]Fig. 12 One-dimensional schematic diagram of the characteristic function in grid[110]

式中：H为特征相函数；x为场中某点的位置矢量；t为时间； ϕ为水平集函数； β为表达插值函数跳跃陡度的经验系数。

采用VOF 法，利用体积组分来计算界面曲率（用于表面张力的计算）比较困难，而采用Level-set方法求解关于符号距离函数的输运方程时，因符号距离函数在界面处连续光滑分布，使得一些关键界面几何信息（最小距离、法线方向、曲率等）的获取更为准确，不过传统的Level-set 方法较难保持流体质量守恒。近年来，CLSVOF 法[111]得到了广泛应用，其很好地结合了Level-set 方法中可以较为准确地捕捉界面曲率以及VOF 方法可以保证质量守恒的优点。尽管这在一定程度上增加了计算成本，但基于该方法求得的相界面曲率光滑性较好，便于实施表面张力平衡模型。

2.1.3 大密度比两相流求解方法

船海结构物水气泡混合流是由空气和水组成的两相体系，二者的密度比值较大，因而给精细化数值模拟带来了一定的困难，特别是对于水气泡混合流中的射流、飞溅等复杂界面变形，或是因强剪切作用而形成的两相界面流动现象，在数值迭代过程中，在气相中有可能会出现较大的伪速度，从而导致界面发生不真实的变形，致使空气速度持续增长，加剧两相流数值求解的不稳定性，甚至是引起迭代发散[112]。Ghods 等[113]指出，动量方程和VOF 方程离散方法的不同会导致动量和体积组分的输运存在相位差，由此引起的界面位置误差会导致在界面附近产生错误的动量输运，而气相因密度较小，容易形成较大的伪速度从而破坏数值计算的稳定性。采用守恒形式离散求解N-S 方程可在一定程度上抑制伪动量的产生。

基于该问题，有学者提出了不同方法用来提高大密度比两相流动模拟的稳定性。Rudman[114]提出了一种从VOF 标量方程中构造通量密度的方法，其使用更加致密的网格来更新体积组分。首先，针对网格上VOF 相函数的分布曲线进行积分，计算得到了控制体边界的相分数值；然后，通过插值得到了边界的密度通量，以用于动量的输运。该方法保证了动量和体积组分以相同的离散方法推进，实现了大密度比二维气泡动力学行为的模拟。Bussmann 等[115]借鉴Rudman[114]的思想，将类似的算法移植到非结构化同位网格中，避免使用更致密的网格构造界面通量，降低了计算量，且应用于二维溃坝流动、大密度比液滴定速传播等问题的模拟也验证了算法的高鲁棒性。Vaudor 等[116]修改了Rudman[114]的原始算法，其基于单网格系统建立了一种通过半网格质量守恒计算质量（密度）通量并用其推进动量的方法，如图13所示。图中，i为网格点序号，Δx和Δy分别为网格在x和y方向的尺度，u为流体速度，ρ 为流体密度，下标f 和c 分别代表面心值与体心值。可在每一时间层面的求解中通过引入密度输运方程来更新密度，采用相同的质量通量推进动量，并在由密度输运方程求解得到的伪密度场获得下一时间步的速度场。因在动量及质量（密度）推进过程中使用的质量通量由相同的数值格式构造，故在每一时间步中能够实现动量及质量的一致推进，则可避免因自由面位置误差而引起的界面附近伪动量或伪速度的产生，提高了数值求解的稳定性。这种动量−质量一致推进的思想近期已被一些学者成功应用于大密度比液滴高速砰击湿润壁面[117]及液体金属涂层过程[118]的模拟中，代表了大密度两相流模拟方法的新进展和突破。

图13 动量−质量一致推进法中的速度和伪密度控制体Fig. 13 Control volume of velocity and pseudo-density in the consistent mass and momentum flux computation method

2.1.4 表面张力模型

船海结构物水气泡混合流中包含有许多微观流动现象，而表面张力对这些微观流动现象起着主导作用，对其中流动不稳定现象的发生和发展具有重要影响。水气界面在表面张力的作用下会收缩并断裂，从而形成小气泡或小液滴。CSF 模型[119]是计算流体力学中常用的表面张力模型，其原理是将表面张力转换为体积力然后加入N-S 方程中进行考虑。表面张力源项的计算式为

式中： σ为表面张力系数； α为流体体积分数； κ为界面曲率，其计算式为

式中，ni为相界面单位法向量。

传统的表面张力模型在模拟流体的表面张力时，流体表面受到的表面张力都是法向力，很难保证流体表面的光滑性。汪欢欢等[120]在CSF 模型的基础上对流体表面施加了切向力，以使流体表面更加光滑，结果显示模拟的流体表面张力效果更好。

相比CSF 模型，尖锐界面力模型（SSF）[121]能够在网格分辨率相对较低的情况下正确模拟该表面的张力效应。将SSF 模型与CLSVOF 方法及虚拟流体法（ghost fluid method, GFM）[122]配合使用，并将表面张力施加在两相界面上，对于毛细不稳定性其有着较高的解析精度，因此可以实现对表面张力占主导作用的小尺度界面流动问题的精细化模拟。可采用CLSVOF 方法计算得到的符号距离函数 ϕ直接计算界面曲率κ：

其中，符号距离函数 ϕ经VOF 函数场重构[123]获得：

式中：Ng为周围交界面网格的数目；αg为网格g中的体积分数值；Dgij为网格g中界面线段至网格（i，j）的法向距离； θgij为网格g中界面线段中心至网格（i，j）中心连接线与界面法向的夹角。

完成曲率计算后，将表面张力以压力梯度算子的形式加入N-S 方程右端，并通过GFM 法引入插值模板点，然后通过界面间断关系构造面心压力梯度，如图14 所示。图中：x为网格空间坐标；x=xin为相界面位置；x¯为构造的虚拟网格点坐标；p为压力；p1(x)，p2(x) 分别为气相和液相的压力。采用该方法，可以提高压力泊松方程及速度修正方程中压力梯度项的离散精度。

图14 虚拟流体法构造的界面压力梯度Fig. 14 Construction of the interface pressure gradient using ghost fluid method

2.1.5 浸入边界法

船海结构物水气泡混合流要考虑壁面边界的作用。在两相流精细化模拟当中，要计算固−液−气三相相互作用的力学特性并模拟交界面的形态比较复杂，且对于船体等复杂曲面外形边界的描述也比较困难，尤其是在考虑航行船舶等运动物体时，随着边界的运动，每一个时间步网格都需要进行更新以适应边界的变化。为了解决这一问题，Peskin[124]提出使用固定笛卡尔网格的浸入边界法，其基本思想是将复杂结构的边界条件模型化为简 N-S 方程中的体积力源项，使用简单的直角网格系统不仅可以回避贴体网格生成困难的问题，还能极大地简化网格生成工作，易于处理动边界问题。

按照对体积力源项的处理方式，可以将浸入边界法分为2 大类，即连续力法（continuous forcing approach）和离散力法（discrete forcing approach）。连续力法最开始主要用于处理弹性边界问题，此时，体积力源项的表达与物面的材料性质有关，由力学关系可以直接得到其表达式。典型的代表是由Peskin[124]提出的浸入边界法。随后，Goldstein等[125]和Saiki 等[126]对该浸入边界法进行了改进，其利用反馈力来表示固体边界，通过指定流体速度在流场某一期望点为0 来满足无滑移边界条件。离散力法的体积力源项由离散后的控制方程计算得到，依赖于控制方程的离散形式，一般无法获得其解析表达式，主要用于处理固体界面问题。典型的离散力法有直接作用力法[127]、切割网格法（cut-cell）[128]以及虚拟网格法（ghost cell）[129]。

浸入边界法提出至今，已从计算精度、时间步长限制、高雷诺数流动等几个方面进行了研究改进。Liu 和Hu[130]提出了一种适用于复杂运动物体周围流场求解的浸入边界法，即运用Level-set方法追踪移动物面：

式中：点P和点A如图15 中所示，其中点P为流场中某点，点A为距离点P最近的物面点，PA为两点间的方向矢量；n为点P附近物面三角形面单元的外法线矢量，图中n1和n2为n的示意举例，nA为A点处的法向矢量。在该研究中，改进了速度边界插值方式，用以抑制压力振荡，并允许在数值求解过程中使用较大的时间步长。随后，又在前期工作的基础上引入了基于移动最小二乘法的高精度边界插值算法[130]，通过结合物体六自由度运动求解技术及网格自适应加密技术，实现了不可压缩流体−刚体耦合运动的精确预报及流场求解。Zhang 等[131]将浸入边界法与几何重构VOF 法（VOF/PLIC）相结合，模拟了二维水平圆柱的入水过程，随后，又进一步将二阶精度的尖锐界面浸入边界法应用于运动固体的模拟[132]中，研究了圆柱入水过程中自由面的变化。Yang和Stern[133]在考虑物面润湿性影响的基础上，将一种简单的接触角边界条件与尖锐界面浸入边界法相结合，模拟二维圆柱入水过程，结果显示，引入接触角边界条件后，圆柱物面后方非物理的闭合空腔消失了，自由面形态更符合真实情况。

图15 用Level-set 方法判断的网格与物面相对位置[134]Fig. 15 Determination of the relative position of grid points and object surface using Level-set method[134]

2.1.6 无网格粒子法

以上介绍的方法均为针对水气泡混合流精细化模拟的网格类方法，虽然采用了诸多技术来克服数值耗散和边界处理等问题，但受限于欧拉求解框架的限制，在复杂运动边界和大变形破碎自由面问题上仍有网格处理上的困难。与之相比，无网格粒子法是基于拉格朗日框架构建的求解算法，该方法将连续流体视为带有速度、压力等物理特性的流体粒子，基于梯度算子和拉普拉斯算子离散控制方程，可以有效避免对流项网格离散求解的耗散问题。同时，无网格粒子法可以通过判断自由表面的粒子来方便地描述自由面的演化，避免了通过网格来描述自由面时的困难[38]。因此，近年来有学者开始采用无网格粒子法开展水气泡混合流相关问题的模拟探索。

文潇[135]发展了一种高精度和高稳定性的多相流移动粒子半隐式（parallelization of moving particle semi-implicit, MPS）方法，其采用各相同时求解的耦合策略，通过引入密度光滑格式、相互作用黏度、表面张力模型以及两相粒子碰撞模型，实现了对复杂多相流动的模拟。然后，在此基础上又面向大密度比的水气泡流动问题，加入了基于粒子间平均密度的压力泊松方程求解、适用于高密度比的压力梯度模型，以及可压缩−不可压缩耦合模型，使得采用该方法能够对水气高度掺混时的界面压力和速度实现准确求解。Wen 等[34]和陈翔等[38]基于该无网格粒子法，对高度分散的微尺度气泡液滴生成演化过程进行了数值模拟，在模拟过程中，采用流体粒子可以很好地模拟相界面翻卷破碎、水跃、射流飞溅、气泡群扩散等强非线性现象，且研究还进一步对因多相掺混引起的流场速度和压力脉动进行了分析。为了加速计算，无网格粒子法可与GPU 计算技术相结合[36]，以克服因流体粒子数目过多带来的计算效率下降的问题。

采用无网格粒子法模拟船海结构物水气泡混合流问题的局限性主要包括2 个方面[136]。第1 个方面是计算效率的问题。对于三维物理问题的模拟，邻居粒子数目的增加会明显增加相互作用搜索和模型处理的时间，同时也会增加压力泊松方程求解的维度与困难程度，致使计算效率明显降低，导致该方法在三维大尺度问题上的应用出现困难。第2 个方面是边界层黏性流动的求解问题。无网格粒子法很难像网格类方法一样布置出致密的边界层网格，采用该方法模拟高雷诺数黏性流动效应和流动分离效应存在一定的缺陷，其对实际的复杂船舶和海洋平台周围流场的模拟存在局限性。

2.2 水气泡混合流模式化模拟方法

船海结构物水气泡混合流问题涉及超大时空尺度跨越，其中微尺度气泡液滴为毫米量级，聚并破碎等动力学行为只需要几毫秒，而宏观尺度水气泡混合流尾迹却可达千米量级，完全生成或是消除需要几十分钟[137]。同时，水气泡混合流广泛存在于实尺度结构物周围，而在模型尺度结构物实验中却不明显，因而有必要开展大尺度船舶海试和模拟研究。基于以上两方面的原因，要针对实际问题进行完全精细化的求解非常困难，因此有必要对其中的大规模微气泡群等小尺度结构进行模式化建模，先通过气泡卷吸模型等描述其生成特性，然后再通过宏观界面捕捉方法与分散两相流模拟方法的结合形成跨尺度耦合方法，以对船海结构物水气泡混合流进行模拟。下面，将从3 个方面介绍其数值方法的研究进展。

2.2.1 气泡卷吸物理模型

气泡卷吸是气泡混合流问题中的一个重要现象，寻求正确的模式化建模方法来模拟气泡卷吸现象，对于研究上述问题十分关键，同时，也是进一步研究气泡混合流与结构物相互作用的基础。

诱导气泡卷吸现象的物理过程主要发生在自由面附近，根据Castro 等[84]和Li[138]的分析，该过程可以分为直接卷吸和涡卷吸2 种。其中直接卷吸包括因波浪破碎而产生的射流冲击水面卷吸气泡入水、由波浪翻卷形成的气穴卷吸气泡入水以及由飞溅液滴撞击水面捕获气泡导致的气泡卷吸等情况，分别如图16（a）～图16（c）所示。对于涡卷吸，因湍流以及其他因素导致的旋涡会使得自由面发生变形，当变形过大时，可能会诱导气泡卷吸现象，气泡此时受到湍流旋涡作用力、浮力和表面张力的联合作用。通常，切向旋涡和垂向旋涡也会诱导不同的气泡卷吸运动，分别如图16（d）和图16（e）所示，相比切向旋涡，垂向旋涡通常卷吸气泡的入水深度更深。图中，ω 为漩涡旋转的角速度矢量。

图16 诱导产生气泡卷吸的物理现象[138]Fig. 16 Physical phenomena that induce bubble entrainment[138]

近年来，学者们提出了许多气泡卷吸模型。Baldy[139]基于相似性原理，对自由面处因湍流诱导波浪破碎产生的气泡卷吸现象进行研究，提出了相应的气泡卷吸模型。该模型基于Clift 等[140]的研究，假设气泡卷吸速率与湍流耗散率ε呈线性关系，研究得到在自由面附近，气泡的聚集程度正比气泡直径d的负二次方，该卷吸模型表达式为：

式中：S为卷吸量； ε¯为单位质量流体的湍流耗散率； ρ为密度； σ为表面张力系数；K1为模型系数。

Shi 等[15]对由波浪翻卷破碎诱导的气泡卷吸进行了研究。其结合Deane 和 Stokes[12]在实验中得到的气泡尺寸分布谱，将自由面处的气泡聚集程度与剪切作用建立了联系。气泡卷吸源项模型可以写为：

式中： εn,i为在气泡数密度对流扩散方程中添加的源项； αb为气体体积分数；Pr=µt|S|2，为剪切生成项，其中 µt为涡黏系数；Di为与气泡直径相关的概率函数；drb,i为气泡半径空间分布。

Shi 等[15]将该气泡卷吸模型应用到了双流体方法当中，结果显示，该模型在参数选取合适的条件下可以捕捉到气泡卷吸过程的主要特点，但对于波浪刚开始破碎时，存在较大气穴，模型未能表现出好的模拟效果。

Ma 等[141]基于Baldy[139]的研究，同样假设气泡卷吸速率与湍流耗散率ε呈线性关系，基于Carrica等[142]提出的多分散双流体方法，建立了气泡卷吸模型，该模型可以写为：

式中：S(ri)为单位体积产生气泡的速率；f(ri)为气泡群粒径谱；Cb为控制卷吸气泡量的系数； ρl为液体密度； αl为液体体积分数；ri为第i组气泡的半径；NG为划分的气泡种群数目。

Ma 等[141]采用求解器TRUCHAS，对波浪冲击海滩翻卷破碎时的气泡卷吸现象进行了模拟，发现气泡尺寸的分布规律与Deane 和 Stokes[12]的实验结果吻合较好，研究显示，直径小于1 mm 的气泡是由射流与自由面的相互作用产生的，而直径大于1 mm 的气泡则主要是由气穴破碎所导致，并且气泡的存在会在一定程度上抑制由波浪破碎导致的湍流。

式（9）所示的气泡卷吸模型在Derakhti 和Kirby[16]的研究中得到了进一步应用。Derakhti 和Kirby[16]在Carrica 等[142]和Ma 等[141]研究的基础上，将VOF 方法和双流体方法相结合，针对深水孤立波波浪破碎过程开展了数值模拟。其采用VOF方法捕捉自由面，双流体方法模拟分散气泡，湍流模型则使用 LES 方法，并考虑气泡对能量耗散和水相与空气相间动量输运的影响，模拟得到波浪破碎时气泡掺混结果与试验结果[143]的对比如图17 所示。图中，x*，z*分别为根据波长无量纲化后在波浪传播方向和水深方向上的距离，自上而下各图分别表示在参考时间（t-tb）为 0.05，0.15，0.30，0.55，0.75，0.95，1.15，1.45 时的结果。从图中可以看到，水气泡混合流的时空分布得到了较好的复现。研究发现，在波浪破碎区域，气泡诱导的能量耗散占总能量耗散的50%以上。

图17 深水孤立波破碎数值模拟与试验对比[143]Fig. 17 Comparison between numerical results and experimental results of deep-water solitary wave breaking[143]

Ma 等[144]进一步修改了该模型，其假设气泡卷吸速率与流体局部湍动能和流体速度梯度的乘积成正比，提出的亚格子气泡卷吸模型如下：

式中：S(x)为单位时间单位体积内引入的气泡总数密度；g为 由重力引起的加速度值；cent为该模型系数； ϕent为界面厚度；un为流体速度在交界面法向处的分量；k(x)为某位置处的局部湍动能。Ma 等[144]采用双流体方法对俯冲射流和水跃这2 个涉及自由面的气泡卷吸现象进行模拟，得到了较好的模拟结果。

目前，大多数气泡卷吸模型都是基于以前的经验和实验结果，对初始时刻和物理边界处的含气率与气泡尺寸分布给予了一定的假设，并将气泡卷吸速率和湍流进行了联系以进行进一步的分析研究。对于不存在结构物壁面而只考虑自由面的气泡卷吸问题时，例如波浪破碎，这些模型体现出了较好的模拟效果；但当考虑壁面存在时，由于结构物壁面诱导的湍流往往远离自由面，其诱导的气泡卷吸问题有所不同，因而以上模型将不适用[84]。

船舶结构物在海上航行时会卷吸气泡，进而形成水气泡混合流。由于有关船海结构物的气泡卷吸问题涉及的尺度较大，采用直接模拟的方法需要较多的计算资源，因此有必要寻求合适的气泡卷吸模型。目前，多是寻求能较好地反映气泡卷吸物理特性的亚格子模型来对气泡卷吸问题进行数值模拟。Moraga 等[145]根据流体局部速度是否超过阈值以及距离交界面的距离的准则，建立了可以确定气泡卷吸发生位置的亚格子气泡卷吸模型，其表达式为：

Moraga 等[145]采用多分散双流体方法求解气泡输运，在一定程度上较好地模拟了DTMB 5 415船因艏波浪破碎导致的气泡卷吸和“雅典娜”号船近船体尾迹中的气泡卷吸现象。Ma 等[72]将文献[144]建立的气泡卷吸模型用于模拟“雅典娜”号直航和定常回转时周围的气泡卷吸现象，并将结果与实验结果进行了对比，结果表明采用该模型可较好地模拟船舶周围流场中空气体积分数的分布情况。Castro 等[84]针对旋涡与自由面的相互作用过程进行建模，建立了气泡卷吸模型，并应用REX 流体力学求解器，采用多分散双流体方法研究了实尺度下“雅典娜”号航行时船体周围的气泡卷吸现象，自由面采用单相Level-set 方法进行确定，其模型框架形式为：

该模型克服了之前模型中气泡卷吸仅由湍流参数决定而与气泡相对自由面距离无关的缺点，并充分考虑了气泡卷吸过程，包括由单涡引起的气泡卷吸、复数涡叠加作用、卷吸气泡群尺寸和水下扩散演化过程等，这也进一步反映了气泡卷吸的力学机理。通过与Tavakolinejad[146]的实验结果进行对比，发现该模型可以较好地捕捉气泡卷吸位置和卷吸程度。Li 等[147]基于Castro 等[84]建立的气泡卷吸模型，采用多分散双流体方法，结合重叠网格技术和Level-set 方法，研究了“雅典娜”号船在转舵过程中流场和旋涡对气泡卷吸的影响，以及气泡输运和尺寸分布特性等问题，得到了较好的模拟效果。

2.2.2 欧拉−拉格朗日跨尺度方法

采用欧拉−拉格朗日跨尺度耦合方法模拟水气泡混合流的基本思想是：在欧拉框架下，采用界面捕捉方法模拟结构物周围兴波等自由面演化现象，在拉格朗日框架下，采用离子追踪方法模拟微尺度气泡的运动，忽略气泡表面的变形。气泡运动求解的控制方程为动力学方程，一般具有如式（13）所示的形式，其中待求解物理量为气泡运动速度，方程右端为各种受力模型。这种耦合方法由于气泡的描述不依赖网格分辨率，因此可以用于模拟亚网格尺度的气泡群演化，相比精细化求解，计算效率明显提高，并且离散相的求解不存在数值耗散，适于计算水下大规模气泡群的体积、运动和分布。

式中：m为气泡质量；v为气泡运动速度；ω 为气泡旋转角速度；u为液体速度；ρl，ρg分别为液体密度和气体的密度；FD，FL，FF，FB，FAM分别为阻力、升力、流体加速度力、浮力和虚拟质量力；CD，CL，CVM分别为阻力系数、升力系数和虚拟质量力系数。

目前，针对该耦合方法的研究大多处于算法的设计和优化阶段，主要应用在一些基本流动的科学问题上，例如射流雾化、气泡羽流等。Herrmann[148]将Level-set 方法与拉格朗日追踪方法相结合构建了耦合方法，其中Level-set 方法可以用于模拟占据足够多网格的连续相界面，通过基于尺度判别的转化算法，流场中不能通过网格合理求解的小尺度两相结构将自动转化为拉格朗日粒子进行模拟。利用该方法，明显提高了射流雾化过程的模拟效果[149]。Patkar 等[150]的耦合方法结合了独特的转化算法与拉格朗日撒点模型，并利用该耦合方法对许多科学问题进行了模拟，包括上升气泡破碎、气泡羽流、空化等，典型的计算结果如图18所示。其中，（g）图通过染色对流场中采用不同方法描述的气泡予以了区分，从中可以看出，采用耦合方法可以有效解决因网格分辨率不足而引起的小尺度气泡不能捕捉的问题。在欧拉−拉格朗日跨尺度耦合方法中，掺混的小尺度气泡、液滴量一般不需要额外进行人工模型干预，但在Hsiao 等[151]的研究中，其通过精细的试验测量和理论分析，建立了亚格子尺度掺混模型，并将其应用到了水平射流问题的耦合方法模拟中，结果显示模拟结果与试验结果吻合较好。

图18 采用Level-set 和拉格朗日追踪耦合方法模拟的多尺度气泡流[150]Fig. 18 Simulation of multi-scale bubbly flow using Level-set/Lagrange coupled method[150]

空化是一类典型的多尺度两相流问题，该现象常出现在船舶螺旋桨和舵等推进系统的作业过程中。Hsiao 等[149-152]在拉格朗日粒子追踪方法中添加了基于Rayleigh-Plesset 方程的空泡求解算法，该方程将气泡半径作为未知量进行求解，可以用于模拟空化过程中气泡的膨胀和溃灭，控制方程如式（14）所示。

式中：R，R0分别为t时刻和初始时刻的气泡的半径；us=uenc−ub, 其中ub为气泡的移动速度，uenc为气泡运动过程中周围环境流体的运动速度；penc为周围环境流体的压力；pg0为初始时刻气泡内的压强；pv为蒸汽压； γ为表面张力系数；µ为流体黏度系数；kc为多变压缩常数。

Hsiao 等[149]基于该方法模拟了水翼和螺旋桨表面的空泡初生，尤其是叶梢部位，典型的模拟效果如图19（图中，P为水翼表面压力）所示。由图可见，水翼表面上的空泡初生为微尺度的小气泡，该现象可采用拉格朗日点来表示；随后，在发展的过程中发生膨胀和聚合，形成大尺度片空泡，该现象可采用界面捕捉方法表示；最后，脱落后再次演化形成扩散的小气泡由拉格朗日点表示。结果显示，利用耦合方法后整个空化过程的模拟效果得到了明显提升。

图19 采用Level-set 和拉格朗日追踪耦合方法模拟的水翼空化流[149]Fig. 19 Simulation of hydrofoil cavitation flow using Level-set/Lagrange coupled method[149]

欧拉−拉格朗日跨尺度耦合方法的关键问题在于能否全面考虑多尺度流动现象之间的转化。早期的一些研究提出的算法是半耦合的[153-154]，即只能考虑拉格朗日粒子在与自由面发生接触时的聚并融合，而由自由面翻卷自然破碎形成的气泡掺混却未得到考虑。这种方法适用于一些流场输入条件为微气泡的流动，如气泡羽流等，对于船海水气泡混合流这种因界面自然演化而形成的多尺度流动则难以模拟。在后续的研究中，又建立了全耦合的欧拉−拉格朗日跨尺度耦合方法。Tomar等[155]通过采用识别计算域中空相分数不为0 的网格，判断了界面模拟至拉格朗日粒子的转化，同时，拉格朗日粒子至界面模拟的转化可以通过二者间的间距来判断。该方法思路简单、直接，为后续研究奠定了基础。随后，又有研究人员在此基础上从转化判据的角度对多尺度耦合方法进行了优化。在Ling 等[156]构建的数值方法中，包含超过43个网格，且与自由面距离超过1 倍直径的气泡或液滴将被转化为拉格朗日粒子。反之，与自由面距离小于1 倍直径的拉格朗日粒子将被转化为界面捕捉方法来求解。在Zuzio 等[157]的研究中，转化气泡所需包含的临界网格数为83个，其采用自适应网格技术，通过结合拉格朗日模式化方法，有效减小了自适应网格的加密等级，扩大了时间步长，提高了计算效率。需要注意的是，以上研究中转化算法的实施需要以气泡尺寸识别为基础，常用的识别算法为泛洪填充算法，涉及到对全场网格的遍历和临近单元的搜索、识别与储存。由于此过程在每个时间步都要执行一次，因而会增加计算时间并给内存带来压力，一定程度上限制了计算效率的进一步提升。Zhang 等[158]针对该问题提出了一种基于曲率的转化算法，通过标模计算，对比说明了该方法在计算精度和效率上的优势。

目前，欧拉−拉格朗日跨尺度耦合方法因兼顾了流动捕捉效果和计算效率的优点，所以在一些小规模的多尺度两相流问题中得到了应用，但在船海工程领域应用得较少。对于实际的船海工程水气泡混合流，微尺度气泡液滴结构的数目在百万、千万级以上，而采用拉格朗日追踪方法，由于存在粒子相互作用碰撞搜索、并行效率等算法方面的问题，导致在离散相数目规模扩大后计算量有明显的增加，也很难用于处理大规模的船海工程水气泡混合流问题，因此仍需进一步采用模式化方法。

2.2.3 欧拉−欧拉跨尺度方法

采用欧拉−欧拉跨尺度耦合方法模拟水气泡混合流问题的基本思想是：对于大尺度的自由面，采用VOF 或Level-set 等界面捕捉方法进行捕捉，而对水下掺混的气泡等小尺度界面，则采用双流体模型进行求解，其中双流体方法在求解过程中将气泡混合流分为了连续相（液体）和离散相（气泡），并采用欧拉方法进行了建模。求解的有关连续相和离散相各自的质量守恒及动量守恒方程为[145-146]：

式中：下标k代表连续相和离散相，取不同的值时代表不同的相；Mk为连续相和离散相之间的动量输运，包括虚拟质量力项、拖曳力项、升力项、压力项和浮力项等；Tk为应力张量， αk为体积分数；ρk为密度；uk为速度；g为重力加速度。欧拉−欧拉跨尺度耦合方法对于大尺度界面和小尺度界面的处理均采用了欧拉观点。

目前，对于该种方法的耦合主要存在2 个方面的难点[159]。首先，界面捕捉方法与双流体方法的求解框架不同。对于界面捕捉方法，两相是通过求解一组动量方程来获得速度场，因此速度场是共享的，但对于双流体方法，速度场是通过求解两相各自的动量方程组得到的，速度场是不共享的。其次，对相分数的理解这2 种方法也不同，对于界面捕捉方法，为了捕捉较为精确的界面，体积分数的扩散效应很小，但是对于双流体方法，体积分数的求解强调了网格的平均效应，扩散效应会更加明显。

欧拉−欧拉跨尺度耦合方法在波浪破碎、液体射流和水跃等不含结构物的水气泡混合流问题的研究中已得到应用。Shi 等[15]基于Buscaglia 等[160]研究的双流体模型和Deane 和Stokes[12]的实验数据，发展了一种VOF 方法与多分散双流体方法相结合的耦合方法，即采用VOF 方法模拟气液界面，采用基于求解气泡数密度的多分散双流体方法模拟气泡群，来对波浪破碎问题进行研究。Yan 和Che[159]建立了关于VOF 方法和双流体方法的耦合模型，如图20 所示，其将所有流体分为了连续液体、具有大尺度界面的大气泡和具有小尺度界面的离散气泡三相，这三相的体积分数分别为 α1， α2， α3，提出采用“体积分数再分配”的方法处理纯流体、大尺度界面以及小尺度界面在网格中共存时产生的矛盾。经研究发现，对于小尺度流动的捕捉，该方法的求解效果与双流体方法一致，而对大尺度界面的捕捉，该方法与VOF 方法的求解效果一致。

图20 多尺度方法中流体相的划分[159]Fig. 20 Division of fluid phases in multi-scale method[159]

对于垂直的液体射流撞击自由面和水平水跃的问题，是不考虑结构物的典型水气泡混合流的物理问题的。Xiang 等[161]建立了VOF 方法与双流体方法的耦合模型，其中双流体方法通过结合多种群方法求解气泡尺寸，研究了水跃问题过程中的气泡卷吸以及气泡尺寸分布的变化情况，其与实验数据[162]的对比显示吻合较好。Ma 等[144]结合双流体方法与Level-set 界面捕捉方法，提出了气泡卷吸模型，并针对垂直液体射流撞击自由面和水跃的问题予以了研究，如图21 所示。图中，α 代表空气的体积分数，X为水流流向方向，Y为水深方向，X和Y均通过上游水深Dj进行了无量纲化。模拟结果与实验结果的对比显示，建立的模型可以较好地模拟相分数分布。

图21 欧拉−欧拉跨尺度耦合方法模拟不含结构物的水气泡混合流问题[144]Fig. 21 Simulation of water-air-bubble mixed flow without structure using Euler-Euler multi-scale coupling method[144]

在欧拉−欧拉跨尺度方法当中，对掺混气泡群的另一种模式化求解思路是以气泡数密度为基础构建控制方程。基于玻尔兹曼理论得到的离散相气泡数密度方程形式如下[142,163-164]：

式中：下标g代 表不同的气泡种群；ug为气泡速度；Ng为气泡数密度，可用于进一步积分获得气泡体积分数； υt为无量纲涡黏系数；S cb为气泡的施密特数；Dp为仅与体积分数有关的扩散系数； βg， χg，Sg分别为气泡因破碎、聚并和卷吸而导致的生成源项。

由下式，可以获得气泡质量的输运，同时，还需要配合描述气泡速度的方程来表达动量守恒。动量方程的形式如下[142,163-164]：

在此基础上结合基于Rayleigh-Plesset 方程的质量输运模型，还可以实现空化过程的模拟[164]。

基于该思想，Carrica 等[142]首先建立了用于研究船舶水气泡流动问题的欧拉−欧拉跨尺度耦合方法求解框架，其采用多分散双流体方法并结合Level-set 方法，对波浪中水面船舶周围不可压缩水气泡流动问题予以了研究。Castro 等[84]在Carria等[142]研究的基础上，发展引入了气泡卷吸模型，考虑到螺旋桨、舵以及其他附体等的存在，结合动态重叠网格技术，研究了“雅典娜”号船在直航情况所产生的气泡混合流问题。Li 等[163]提出了一种可以改善求解气泡混合流问题的鲁棒性方法，改善了压力速度耦合算法，可以实现高含气率下的水气泡混合流问题，成功对实尺度下一艘平底船周围的水气泡混合流现象予以了模拟。Li 和Carrica[165]考虑到气泡的可压缩性，提出了可压缩气泡混合流问题的求解方法，并模拟了“雅典娜”号船在2 个不同航速下的气泡卷吸现象，结果发现船舶航行中，在船体表面、船艉以及船艏波浪破碎处会有较多的气泡混合流现象发生。相对于高速情况，在低速情况下艉部的气泡卷吸将更加明显，含气率也会更高。Li 等[147]基于以上研究方法，将其应用到了实尺度“雅典娜”号船在Z 形操纵运动下的水气泡卷吸现象研究中，如图22（图中，Sg为气泡卷吸源项大小）所示，得到了其在操纵过程中船体周围水气泡混合流的密度、含气率、涡结构等特征，研究发现在操纵过程中，船艉处的气泡卷吸程度会加剧，但在船体表面会减弱，其尾迹相对直航情况下也会更宽。

图22 欧拉−欧拉跨尺度耦合方法模拟“雅典娜”号船在操纵运动下的艏艉含气率分布[147]Fig. 22 Simulation of bow and stern void fraction distribution of Athena ship in maneuvering using Euler-Euler multi-scale coupling method[147]

相对于直接数值模拟方法和欧拉−拉格朗日跨尺度耦合方法，由于引入了更加模式化的方法求解气泡群，欧拉−欧拉跨尺度耦合方法属于多相流问题中的宏观建模方法，在求解过程中不对气泡的具体界面和形状进行捕捉，而是着眼于捕捉流场中离散相的相分数分布。因此，此方法的求解效率会更高，对于求解实际船海结构物水气泡混合流这种跨尺度问题的研究具有潜在优势。

3 总结与展望

水气泡混合流是一类在船海结构物周围广泛存在的复杂流动现象，其包含水、气、泡、雾等多种流动状态，具有多相、多尺度、多物理场等特点。目前，主要的研究重点在流动机理和数值模拟方法两方面，本文以这2 个核心问题为主线，总结凝练了当前关键的研究方向和研究进展，目标是为研究人员认识、理解船海结构物水气泡混合流机理提供基础，为进一步形成有效的预报、分析和设计技术提供参考。

在流动机理方面，以从微观到宏观的角度系统介绍了目前船海水气泡混合流机理的研究现状。在微观上，总结了当前以波浪破碎和气泡/液滴等微观界面流动为对象的水气泡混合流生成演化机理探索工作，从局部水体多相界面撞击角度介绍了自由面破碎和气泡、液滴结构的生成机制，进一步介绍了生成后的微观界面演化机理研究。从黏性和表面张力这2 个方面阐述了当前对水气泡混合流动影响因素和影响效果的讨论。最后，介绍了水气泡混合流对噪声源和声波辐射的影响研究。在宏观上，以船体和平台周围不同位置流动差异为线索，分类阐述了水气泡混合流现象机理及其对结构物性能影响的现有结论，包括船艏、船舯、船艉型线导致的水气掺混机理差异和演化区别，介绍了水气泡混合流对平台砰击和上浪性能的影响研究。这有助于研究人员更好地认识和理解该现象，拓展研究思路。

在数值模拟方法方面，从精细化和模式化2 个角度，详细介绍了适用于水气泡混合流研究的数值模拟方法。在精细化方法上，首先针对网格类方法从求解精度、效率、稳定性等几个方面介绍了AMR 方法、质量/动量一致推进、尖锐界面力模型、浸入边界法等技术手段。对于无网格粒子法，介绍了多相流MPS 方法在剧烈自由表面流动和水气泡混合流上的研究应用。在模式化方法上，介绍了当前用于模式化大规模微观结构的物理模型以及基于欧拉−欧拉观点和欧拉−拉格朗日观点的模式化方法。在以上数值方法介绍中，分析了当前算法的优、缺点和能力，为研究人员进一步开展数值方法研究奠定了基础。

随着船海结构物不断朝着高速化、大型化的趋势发展，水气泡混合流的影响问题也日益凸显，在绿色节能、稳定安全、降噪隐身等方面体现出越来越强的研究需求，对未来的研究展望可总结为以下几点：

1） 针对实尺度船海结构物的水气泡混合流高效、可靠模型化预报方法。根据现有的观测研究，船海水气泡混合流的一个重要特点是在实尺度结构物周围广泛存在而在模型尺度却难以出现。在相似规律不明确的现状下，亟需发展针对实尺度船海结构物的水气泡混合流高效、可靠模型化预报方法，建立能反映真实流动效应的物理模型，解决计算量、预报精度等关键问题，为真实流动分析提供手段。

2） 水气泡混合流对船舶快速性能的影响机理。目前，船海结构物水气泡混合流机理研究的重点在于对流动形态和生成演化过程的分析，对船舶具体性能影响的研究较少。当前，已有部分研究针对水气泡混合流生成演化过程中的能量传递进行分析，但其消耗船舶能量以及引起阻力变化的定量化结论还有待阐明。另一方面，船舶尾流和船侧的气泡下扫可能导致气泡群进入伴流区域，需要系统研究气泡下扫量以及其对螺旋桨推进、空化性能的影响机理。

3） 水气泡混合流对结构物砰击上浪性能的影响机理。结构物砰击过程本身由于瞬时性、极端性等特点具有一定的研究难度，而对纯水波浪砰击和含有大气腔的波浪砰击问题则有一定的研究。根据实际观测，水气泡混合流是结构物砰击上浪过程的重要流动形式，需要进一步探明多尺度气泡、液滴群对砰击载荷上浪载荷所产生的非线性影响，研究其引起的砰击上浪范围、幅值的变化特性。

4） 水气泡混合流与噪声的多物理场耦合机理。气泡、液滴的振荡、聚并和破碎行为以及其与自由面的撞击等可被视为重要的噪声源。船海结构物水气泡混合流中包含数目庞大的气泡液滴结构，对水面舰船流动噪声的影响是未来的一个重要研究方向，需要构建两相噪声预报技术，研究水气泡混合流产生的船舶近场噪声源分布，并在此基础上分析远场声辐射特性，掌握水气泡混合流与噪声的耦合机理。