许海雨，罗 凯，黄 闯，左振浩，古鉴霄

(西北工业大学 航海学院，西安 710072)

超空泡减阻技术通过生成包裹航行器的超空泡，避免了大部分航行体表面直接沾湿，从而实现航行阻力的显著减小和航行速度的大幅提高[1].通气超空化问题受到国内外研究学者的高度重视.通气超空泡通常分为两种泄气模式：回射流泄气和涡管泄气.前者通常发生在弗劳德数较大、重力因素可忽略的情况下[2]，后者通常出现在弗劳德数较低，重力对超空泡流型影响较大的情况下[3-4].众多学者研究得到了双涡管泄气超空泡与回射流泄气超空泡的转变条件[5].

Skidmore[6]在试验中发现航行器模型在小尺寸水洞中空泡以双涡管模式泄气，而在大尺寸水洞中空泡以回注射流模式泄气.Karn等[7]通过试验研究了不稳定来流对超空泡形态及泄气方式的影响，指出闭合区泡内外压差是引起空泡泄气方式转变的原因.Rashidi等[8]通过水洞试验及数值仿真方法研究了通气超空泡的回射流现象，捕捉到空泡内气体向空化器方向流动，并发现空泡内有明显气体旋涡.王志等[9]采用离子图像测速(PIV)技术对通气超空泡内气体流动结构开展研究，发现了3种典型的空泡流动结构并对每种结构的成因进行了分析.Wang等[10]采用试验与数值模拟相结合的方法建立了通气超空化数值模型，模拟了涡管泄气和回射流泄气超空泡，并对其流动结构开展了研究.张孝石等[11]在水洞中开展了通气量对航行体云状空泡稳定性影响的试验，指出通气量影响空泡演变过程，通气量增大导致空泡发展周期增大，空泡脱落尺寸减小.Ahn等[12]在两种不同水洞设备中对比分析了相同工况下的超空泡外形及尺度的变化规律.Cao等[13]采用欧拉多相流模型研究了双涡管泄气超空泡的空化特性，得到了空泡形态及泡内压力随通气率、流速的变化规律.邓飞等[14]在水洞试验中研究了双圆盘空化器的超空泡形态变化特性，指出双圆盘空化器存在前盘优先诱导空泡生成及后盘诱导空泡生成两种空化状态.Xiang等[15]对通气超空化流动特此开展数值仿真研究，得到了空泡闭合区气体体积分数、气泡速度及气泡尺寸的分布特性.

综上，上述研究均集中于发展充分的双涡管泄气或回射流泄气超空泡形态及流场特性，缺少对超空泡初生及发展过程中空泡演变特性的研究，特别是关于在空泡发展过程中空泡形态及泄气方式的动态变化特性研究.本文基于分相流模型和SST湍流模型，建立可压缩通气超空化数值模型.在通气率和弗劳德数(Fr)系列变化范围内，系统地研究两种超空泡的生成及发展演变过程，重点分析了空泡生成过程中空泡形态及泄气模式的变化规律.

1 计算方法

1.1 数值方法

1.1.1控制方程 商用流体计算软件CFX的多相流模型包含匀质平衡流模型和分相流模型，文献[13]分析了两种多相流模型对通气超空化的影响，指出分相流模型在计算通气超空化上具有较高的计算精度.

通气介质为常温空气，忽略两相间换热，因此能量方程可以忽略，控制方程只涉及连续性方程、动量方程及体积分数方程.

连续性方程为

(1)

式中：γα为α相(1为水相，2为气相)的体积分数；ρα为α相的密度；ui为速度；i和下文j分别表示坐标方向；t为时间.

动量方程为

(2)

式中：p为静压；μ为动力黏度；Mi为相间作用力.

1.1.2湍流模型 基于Baseline (BSL)k-ω模型的SST湍流模型考虑了湍流剪切应力的传输，可以精确地预测流体流动的开始和大逆压梯度条件下的流体分离.因此本文采用SST湍流模型，基本方程为

(3)

(4)

式中：ρ为密度；k为湍流动能；μ为流体动力黏度；μt为涡黏度系数；pk为湍流生成率；ω为湍流频率；σω2、σω3、σk3、α3、β3和β′为模型常数.

(5)

(6)

(7)

式中：ν为运动黏度；y′为到最近壁面的距离.

1.1.3算法配置 本文采用有限体积法离散控制方程，控制方程的瞬态项及扩散项分别采用二阶向后欧拉式和中心差分格式进行离散，此外，控制方程的对流项采用高精度格式离散方法以提高计算精度.采用全隐式多网格耦合求解技术，同时求解连续性方程、动量方程及能量方程等控制方程，可进一步提高计算收敛速度.为了节约时间，获得合理空泡演变规律，对非定常计算的时间步长进行无关性分析，最终确定时间步长为0.1 ms.

1.2 物理计算模型

1.2.1计算模型 为了在水洞试验中验证数值方法的准确性，本文采用两套计算模型，一套为水洞试验验证的细直杆模型；一套为研究通气空泡生成、发展的数值计算模型.为了减小模型对通气超空泡流型的影响，水洞试验采用细直杆模型，如图1(a)所示.数值计算模型为了清晰捕获空泡生成及发展过程中的空泡形态演变特性，减少模型对其影响，采用单空化器模型，如图1(b)所示.

图1 几何模型Fig.1 Geometric model

1.2.2计算域及监测位置设置 通气超空化区伴随空泡脱落及强湍流度，因此需要划分足够精细的计算网格，对整个计算域进行结构化网格划分，其中60%以上的网格分布在两相空化区.并对计算网格进行网格无关性检验，网格数为80万.计算域左侧边界采用速度入口，计算域右侧边界采用压力出口，通气孔采用质量流量入口边界条件，通气为可压缩常温空气.计算域周围采用无滑移壁面边界条件.为了研究超空泡内压力的变化特性，分布选取P1和P2作为监测点，以监测流场压力时变特性，其距离空化器的位置分别为2dc和7.5dc，其中dc为空化器直径.如图2所示，图中x为超空泡发展方向，y为重力反方向，来流速度分3种工况：5 m/s、7 m/s及9 m/s.

图2 流域及边界条件Fig.2 Flow domain and boundary conditions

1.3 模型验证

为了评估仿真结果的置信度，在西北工业大学水洞中开展用于数值计算方法验证的水洞试验，数值仿真采用与试验相同的计算工况.图3给出了两种通气率下数值计算与水洞试验的空泡形态结果对比，图中：Lmax为空泡最大长度；Dmax为超空泡最大直径；Cq为通气率.两种工况的弗劳德数相同，均为24.8，Cq定义如下：

(8)

式中：Q为超空泡气体体积流量；v为来流速度.

可以看出，在小通气率工况下，数值计算得到的超空泡最大直径与水洞试验结果误差约为2.82%，超空泡最大长度与水洞试验结果误差为2.75%.大通气率工况下，数值计算得到的超空泡最大直径与水洞试验结果误差约为2.65%，超空泡最大长度与水洞试验结果误差为2.66%.为了进一步验证所建立数值模型的计算精度，图4给出了在两种通气率下数值计算得到的不同时刻超空泡最大长度与水洞试验的对比结果.由图可知，小通气率时，超空泡随时间发展迅速，空泡长度增长较快，随后空泡发展缓慢，超空泡长度的增长率显著减低，最后超空泡发展稳定，超空泡最大长度几乎不随时间改变.这是由于开始时超空泡以回射流方式发展，空泡发展迅速，当流动时间为0.9 s时，空泡泄气方式转变为双涡管泄气方式，由于一部分气体从涡管流出空泡体，因此空泡的发展缓慢.数值仿真可以精准模拟空泡的发展演变特性，并且在不同时刻下空泡尺度与试验结果吻合度较高.在大通气率下，本文所建立的数值模型仍可以高精度地模拟通气超空化的演变特性，在不同时刻空泡长度与试验误差的最大值为2.4%.综上所述，本文所采用的数值计算方法可以有效地揭示通气超空化流动机理，空泡尺度与尾部泄气方式与试验结果符合度较高，这说明本文的计算方法有较高计算精度，可用于通气超空泡生成及演变特性的研究.

图3 Cq不同时数值仿真与水洞试验结果比对Fig.3 Comparison of test and numerical results at different Cq values

2 计算结果及分析

空化数σ及弗劳德数Fr两者决定了通气超空泡的泄气类型，而通气率是影响空化数的最主要参数之一.因此本文选用通气率及弗劳德数两个参数变量研究通气超空泡泄气转变特性.

2.1 通气率对空泡初生及发展演变特性的影响

为了研究通气率对超空泡泄气类型的影响，本文在通气率从小到大范围内，选择一系列通气率开展空泡初生及发展演变特性研究，发现随着通气率的增大，超空泡泄气模式依次为回射流泄气、空泡泄气转变临界状态及双涡管泄气.为了对比研究3种泄气类型空泡的初生及发展演变特性，文中分别选择代表3种典型空泡泄气方式的通气率(Cq=0.03，0.06，0.39).

图5 不同Cq下超空泡轮廓线Fig.5 Supercavity profiles at different Cq values

图5给出了超空泡对称面上的空泡轮廓曲线.由图可知，随着通气率的增加，超空泡整体几何尺寸增大，并且重力对空泡形态的影响也越来越显著.通气率为0.03时，空泡形态较小，空泡长度为5.1dc，并且空泡几何形态为轴对称空泡，几乎不受重力效应影响，超空泡以回射流模式泄气.通气率增加到0.06时，空泡尺寸显著增大，空泡长度为10.38dc，超空泡泄气模式为双涡管泄气，但涡管形态较小.当通气率增加到0.39时，空泡长度增加到26.3dc，超空泡流型受重力影响更大，空泡几何中心线向上弯曲更显著，空泡为非对称几何形态，此时超空泡以双涡管模式泄气.

图6给出了通气率为0.03时的超空泡初生、发展及泄气模式演变过程.由图可知，在空泡刚生成(t=0.04，0.08 s)时，空泡发展非常不稳定，出现了空泡气团的断裂、分离现象.随着空泡的进一步发展，空泡断裂现象逐渐消失，但空泡形态具有较强的非定常特性.超空泡泄气模式为回射流泄气，在空泡闭合区通入的非凝结气体裹挟液体水流向空化器方向，并且在闭合区周期性地发生空泡气团的脱落与空泡的再发展，空泡气团脱落主体空泡后空泡尺寸减小，然后再继续发展变大.例如t=1.62 s时刻空泡气团脱落，空泡长度为46dc，当t=1.66 s时，空泡长度增大到55dc.当t=1.76 s时，空泡气团再次脱落，空泡长度减小到47.5dc.

图6 Cq=0.03时超空泡发展过程Fig.6 Development of supercavity at Cq=0.03

图7为空泡内压力振荡特性，图中pc为空泡内压力.在空泡初生时，P1和P2位置的压力先开始逐渐增大，最大压力分别为97.1 kPa，99.3 kPa.这是因为超空泡刚发展时，监测点位于空泡闭合区后，其压力受到超空泡闭合区的高压影响.随着空泡的进一步发展，当监测点位于空化区时，泡内压力开始减小，但由于回射流现象导致空泡内压力振荡，空化数发生改变，因而空泡形态及尺度也相应改变.当空泡发展充分时，P1位置的压力振荡幅度为3.3%，而P2位置压力振荡更加剧烈，振荡最大幅度为13.3%，其原因在于通气率过小，空泡形态较小，P2正位于空泡闭合区附近，因此压力振荡显著.由空泡形态及压力特性可知，回射流超空泡非常不稳定，难以预估其空泡形态及尺寸.

图8为Cq=0.06时通气超空泡生成及发展演变过程.图中将超空泡在重力方向投影所得空泡形态作为正视图，正视图中空泡尾部出现的管状形态即为双涡管泄气模式，这可以从图中最后一行给出的俯视图中看出.通气超空泡在刚生成时空泡形态表现出强烈的非定常特性，空泡以回射流模式泄气.空泡发展非常不稳定，仍出现了空泡气团的断裂、分离现象(t=0.04 s).随着空泡的进一步发展，空泡逐渐增大，在空泡闭合区因回射流强度及回流方向具有时变性，当气/水混合物与空泡轮廓相互接触时，空泡形态被破坏，呈现出形态各异的超空泡流型(t=0.18,0.32,0.48 s).超空泡在发展充分时空泡泄气方式在回射流泄气和涡管泄气两者之间交替进行，回射流泄气空泡形态如图中t=1.62,2.10 s所示.涡管泄气空泡形态如图中t=2.00,2.32 s所示.

图7 Cq=0.03时空泡内压力振荡特性Fig.7 Oscillation characteristics of pressure inside supercavity at Cq=0.03

图8 Cq=0.06时超空泡发展过程Fig.8 Development of supercavity at Cq=0.06

空泡内压力振荡特性如图9所示.由图可知，空泡内压力发生一定程度的振荡，并且P2位置的压力振荡更为剧烈.按照Campbell-Hilborne空泡泄气转变经验判断准则，σFr<1时，空泡以双涡管模式泄气；当σFr>1，空泡以回射流模式泄气.由于空泡内压力持续振荡，空泡内压力处于最大时，计算得到的σFr=0.96，当空泡内压力处于最小时，σFr=1.58.因此空泡泄气模式不稳定，在回射流泄气和涡管泄气之间来回转变.

图9 Cq =0.06时空泡内压力振荡特性Fig.9 Oscillation characteristics of pressure inside supercavity at Cq=0.06

图10 Cq=0.39时超空泡发展过程Fig.10 Development of supercavity at Cq=0.39

图10为Cq=0.39时空泡生成及发展演变过程，由图可知，通气超空泡在刚生成时空泡发展速度较快.相比较于Cq=0.06的工况，在相同时刻时，超空泡尺度更大，发生空泡形态分离、断裂等现象减少，并且空泡发展速度更快，表明增大通气率抑制了空泡的破裂，有助于超空泡的生成及发展.超空泡开始先以回射流泄气方式生成及发展，如图中t=0.18 s之前的空泡形态.随着空泡的继续发展，空泡几何尺寸增大，受重力效应的影响也越来越显著，超空泡泄气模式最后转变为回射流模式泄气(t=0.8，1.2 s).

图11为空泡发展过程中泡内压力振荡曲线，由图可知，开始时空泡内压力振荡剧烈，超空泡非定常特性强烈，这是由于空泡受到回射流作用导致的.当t>0.8 s时，空泡内压力振荡特性减缓，空泡最大压力振荡幅度仅为344 Pa，此时超空泡以双涡管模式泄气，这表明相比较于回射流泄气，双涡管泄气超空泡更加稳定，空泡内压力及空泡形态非定常性较弱.当超空泡以双涡管泄气时，空泡内无回射流现象，通入气体沿着空泡壁向下发展，在空泡闭合区从双涡管流出空泡体.

图11 Cq=0.39时空泡内压力振荡特性Fig.11 Oscillation characteristics of pressure inside supercavity at Cq=0.39

图12 Fr=11.3时超空泡发展过程Fig.12 Development of supercavity at Fr=11.3

2.2 Fr对超空泡生成及发展演变特性影响

本文选用Cq=0.06，Fr=15.8超空泡泄气转变临界状态作为参考点，通过增大和减小Fr，以获得其对超空泡泄气模式转变的影响特性.因此保持通气率0.06不变，系列改变Fr(分别取11.3、20.3)，以研究其对通气超空泡泄气类型的影响特性.图12为Fr=11.3时的空泡演变过程，在空泡开始生成及发展过程中，空泡主要以回射流泄气模式发展，其空泡演变特性与Fr=15.8时的工况相同，但发展过程中出现涡管泄气次数更多(t=1.06，1.52 s).当空泡充分发展时，空泡泄气方式转变为双涡管泄气，如图中最后一行所示.此时双涡管管径较小，但空泡内部仍存在回射流现象，即空泡尾部气/水混合物流向空化器方向.

空泡内压力振荡特性如图13所示.当空泡泄气类型转变为双涡管泄气时，空泡内压力最大振荡幅度为690 Pa，相比于通气率0.39时工况，压力振荡幅度增大，这也表明空泡内回射流现象的存在.通气率0.39时，空泡形态较大，因此受重力效应影响大，泄气方式完全转变为涡管泄气，并且空泡内气体完全从两个涡管泄出.而Fr=11.3工况空泡形态较小，空泡未全部转变为涡管泄气，因此涡管管径较小，泡内存在回射流动.

图13 Fr=11.3时空泡内压力振荡特性Fig.13 Oscillation characteristics of pressure inside supercavity at Fr=11.3

图14为Fr=20.3时的通气超空泡生成、发展及泄气特性图.由图可知，Fr增大时，通气超空泡泄气模式为回射流泄气，超空泡在生成及发展过程中，回射流对空泡轮廓造成的破坏减小，并且空泡发展速度更快.相比于Fr=15.8工况，超空泡涡管泄气模式消失，在整个发展过程中，超空泡始终以回射流方式泄气(t=2.37，2.40，2.50 s).另外，P2位置的压力要高于P1位置的压力，这表明空泡内压力并非恒定的，具有非均匀特性，泡内压力沿空泡发展方向逐渐增大.

图14 Fr=20.3时超空泡发展过程Fig.14 Development of supercavity at Fr=20.3

图15 Fr=20.3时空泡内压力振荡特性Fig.15 Oscillation characteristics of pressure inside supercavity at Fr=20.3

空泡内压力时变特性如图15所示.相比于低速工况，当超空泡发展完全时，Fr增大，空泡内压力振荡特性增强，压力振幅及振荡周期均显著变大，这与空泡回射流泄气及泡内流动结构有关.通气率为0.03时空泡泄气方式为回射流方式，空泡内压力振荡更加猛烈.由此可知，回射流超空泡内压力具有振荡特性，并且增加流速有助于抑制回射流泄气超空泡的压力振荡.

3 结论

本文对回射流泄气及双涡管泄气的空泡生成、发展及空泡形态变化特性进行非稳态数值计算研究，重点对比分析了不同泄气类型空泡的形态变化规律及泄气模式之间相互转变特性，得到以下结论：

(1)回射流泄气超空泡内存在气体裹挟水反向流动的回射流，空泡尾部伴随气团的周期性脱落，在空泡初生时伴随空泡的断裂及分离现象，空泡形态难以预估，泡内压力振荡剧烈.

(2)双涡管泄气超空泡相对比较稳定，空泡形态及泡内压力变化幅度较小，空泡初生时先以回射流模式泄气发展，随着空泡进一步发展，最终空泡泄气模式转变为双涡管泄气.

(3)当空泡处于或接近空泡泄气转变临界时，空泡在初生时以回射流模式泄气，当空泡充分发展后空泡泄气类型一直在涡管泄气和回射流泄气之间相互转变，流场及空泡形态变化特性更加复杂.