CN111709196A - 一种水下高速航行体抗空化空蚀能力评估方法 - Google Patents

Abstract

本发明属于气泡动力学技术领域，具体涉及一种水下高速航行体抗空化空蚀能力评估方法。本发明研究了空化核在水下高速运动的航行体诱导流场中的运动特性，通过数值计算获得流场的速度和压力，揭示了空化空蚀的力学机理。根据航行体的实际几何建立模型，通过观察单个气核在航行体附近的脉动、迁移以及射流等行为，评估此种航行体的抗空化空蚀能力，为设计具有更好抗空化空蚀能力的航行体结构提供指导。

Description

一种水下高速航行体抗空化空蚀能力评估方法

技术领域

本发明属于气泡动力学技术领域，具体涉及一种水下高速航行体抗空化空蚀能力评估方法。

背景技术

水中大量微米级的气核在流过高速航行结构或高速旋转机械低压区附近时，由于环境压力的降低和表面张力作用产生爆发式的生长形成空化气泡群，称为水利空化现象。空化现象不仅降低螺旋桨推进的效率，增大航行体的阻力，并且空化气泡随着流动向下游移动的同时，在非均匀流场和结构边界的作用下将产生非球状脉动，在坍塌溃灭阶段可能产生指向结构表面的高速水射流，气泡射流砰击载荷被认为是造成螺旋桨剥蚀的重要成因之一。一个多世纪以来，人们对静水中壁面附近气泡的脉动特性和射流规律有较为详尽和深入的研究，但是对非均匀流场中迁移的空化气泡的产生机理与运动特性的认识仍有不足。在空化水洞中进行实验是最直接的方法，但是由于相关实验开展难度大，花费高，迄今为止，公开发表的相关实验数据很有限。幸运的是，随着计算机技术的发展，高精度模拟气泡脉动的数值方法出现，为间接研究空化气泡提供了简单廉价的途径。然而，目前非均匀流场、复杂结构耦合作用下单个气泡的运动特性的数值探究仍然很少，因此，有必要开发一种探究水下高速航行体引起的成核空化的数值方法，旨在揭示空化初生和空化演变的机理。

由于实验条件的限制，从实验中得到有限的数据还不足以揭示空泡运动背后的许多力学机理，结合相应的数值仿真能够揭示空泡更多的物理特性，帮助研究人员更好地认识其中的物理规律。现今气泡动力学中的数值研究方法分两大类，第一类是域离散方法，对整个流域进行离散同时采用各种界面追踪/捕捉法获取气液交界面。第二类数值方法是边界离散方法，仅需对流场边界划分网格单元但仍可获得流场内部信息。前者的计算量庞大，距离工程应用还有较远的距离，后者由于问题减少了一个维度，计算效率十分可观。然而，目前边界元方法在游移空泡形成机理的数值研究较少，对于气核在结构绕流流场、结构边界耦合作用下的运动特性核载荷特性的认识仍有不足。因此有必要提出一种预报水下高速航行体引起的成核空化的三维数值方法，该方法可以对流场中空化核运动与载荷进行全面的研究，并与实验相互印证进而提高现有的研究技术和手段。

发明内容

本发明的目的在于提供一种水下高速航行体抗空化空蚀能力评估方法。

本发明的目的通过如下技术方案来实现：包括以下步骤：

步骤1：输入待评估的水下高速航行体，设定气核的初始半径R_initial，对水下高速航行体的几何尺度与气核的初始半径R_initial进行无量纲化处理，取一组特征量对压力、速度、速度势、时间和密度进行无量纲处理，建立具有普遍意义的无量纲化系统和模型；

步骤2：建立已知边界的均匀来流流场，将无量纲处理后的水下高速航行体模型置于该流场中；设定计数器i和最大计算次数，初始化i＝0；

步骤3：获取当前时刻无气核的情况下，均匀来流流场中静止的水下高速航行体表面的压力分布；

步骤4：以水下高速航行体结构剖面中首端位置为原点，以距首端的无量纲轴向距离X为横轴，以距首端的无量纲径向距离R为纵轴，绘制无气核的情况下，水下高速航行体结构剖面绕流表面压力系数C_p曲线；

其中，P_∞为流场中远方来流的压力；U₀为流场中均匀来流的速度；ρ为流场中流体的密度；P_l是水下高速航行体受到的局部表面压力，根据均匀来流中静止的水下高速航行体表面的压力分布获取；

步骤5：获取无气核的情况下，水下高速航行体结构剖面绕流表面压力系数C_p曲线中压力系数C_p最小位置处的坐标(X₀,R₀)，在该位置附近生成初始半径为R_initial的气泡；

步骤6：令i＝i+1，计算流场的扰动速度势φ_d；

其中，p表示流场控制点；q表示边界积分点；S包括所有流场边界；n为边界法向量，指向流场外为其正向；λ为p点观察流场的立体角；G为格林函数；

步骤7：计算气泡溃灭时产生的射流速度u_i和第i次计算中气泡的半径R_i；

u_i＝▽(U₀X₀)+▽φ_d

其中，Δt_i为第i次计算与第i+1次计算之间的时间步长；S_b为气泡的表面积；

步骤8：将计算得到的气泡溃灭时产生的射流速度u_i量纲化，判断气泡溃灭时产生射流的速度和方向；若射流速度大于20m/s，并且方向指向水下高速航行体结构，则判定产生严重的空蚀效应，水下高速航行体的抗空蚀能力不合格；

步骤9：判断是否达到设定的最大计算次数或者已经产生射流抨击，若未达到设定的最大计算次数且未产生射流抨击，则返回步骤6；否则，执行步骤10；

步骤10：获取计算结果中最大的气泡半径；若最大气泡半径小于2.5，则判定水下高速航行体的抗空化能力不合格；若最大气泡半径小于2.5，则判定水下高速航行体的抗空化能力合格，且最大气泡半径越小，抗空化性能越好；

步骤11：输出水下高速航行体的抗空化空蚀能力评估结果。

本发明的有益效果在于：

本发明研究了空化核在水下高速运动的航行体诱导流场中的运动特性，通过数值计算获得流场的速度和压力，揭示了空化空蚀的力学机理。根据航行体的实际几何建立模型，通过观察单个气核在航行体附近的脉动、迁移以及射流等行为，评估此种航行体的抗空化空蚀能力，为设计具有更好抗空化空蚀能力的航行体结构提供指导。

附图说明

图1为高速航行体剖面和表面压力分布示意图。

图2为轴对称与三维空化气泡半径时历变化对比曲线图。

图3(a)为空化核在绕流流场中膨胀数值计算结果图。

图3(b)为空化核在绕流流场中收缩数值计算结果图。

图4(a)为气泡周围的压力场计算结果图。

图4(b)为气泡周围的速度场计算结果图。

图5为预报水下高速航行体引起的成核空化三维数值方法计算流程图。

具体实施方式

下面结合附图对本发明做进一步描述。

本发明适用于预报背景流场中、结构附近的气泡运动及载荷特性，能够预测水中气核在高速航行结构附近局部低压区域成核空化的物理过程，并且在水利机械空化、水下爆炸、微流控、医学、化工等领域具有广泛的应用前景。

本发明研究了空化核在水下高速运动的航行体诱导流场中的运动特性，通过数值计算获得流场的速度和压力，揭示了空化空蚀的力学机理。根据航行体的实际几何建立模型，通过观察单个气核在航行体附近的脉动、迁移以及射流等行为，评估此种航行体的抗空化空蚀能力，为设计具有更好抗空化能力的航行体结构提供指导。气核空化的临界半径为2.5，若计算气泡半径大于2.5证明在该工况下气核将产生爆发式生长形成空化。抗空化能力强的航行体将抑制气核的膨胀，在相同工况下不同的航行体结构形式通过比较气核的最大半径可评价航行体的抗空化性能，最大气泡半径小于2.5且越小则性能越优。气泡溃灭产生射流对表面产生抨击，其速度和方向决定着航行体表面空蚀的产生，如射流速度大于20m/s并且指向结构将对结构产生严重的空蚀效应，通过本发明提出的数值方法可以提取抨击点的速度，并且预测射流抨击的方向。射流速度越小其造成的危害越小，比较射流的速度也可为航行体防空蚀性能设计提供指导。

本发明求解航行体周围流场的压力，可以用于评估航行体表面受到的载荷，为航行体局部防护设计提供指导。航行过程中头部区域往往会产生数倍于大气压的压力，本发明所用航行体轴向坐标为0-0.25的区域内承受最高4倍大气压，需要在此高压区内加强防护。同时本发明提出的数值方法，可以指导空化水洞中进行的单个激光气泡实验，由于实验技术和水平的限制，实验中难以测量的压力等数据可通过本发明的三维数值方法计算，进而揭示实验中发现的新现象背后的机理等。现今气泡与航行体边界层的耦合作用的研究不足，本发明的基于气泡与航行体非边界层流动的耦合作用，可为开展边界层对气泡的耦合作用研究提供对比。因此，本发明具有重要的工程应用价值和学术意义。

一种水下高速航行体抗空化空蚀能力评估方法，包括以下步骤：

步骤1：输入待评估的水下高速航行体，设定气核的初始半径R_initial，对水下高速航行体的几何尺度与气核的初始半径R_initial进行无量纲化处理，取一组特征量对压力、速度、速度势、时间和密度进行无量纲处理，建立具有普遍意义的无量纲化系统和模型；

步骤2：建立已知边界的均匀来流流场，将无量纲处理后的水下高速航行体模型置于该流场中；设定计数器i和最大计算次数，初始化i＝0；

步骤3：获取当前时刻无气核的情况下，均匀来流流场中静止的水下高速航行体表面的压力分布；

步骤4：以水下高速航行体结构剖面中首端位置为原点，以距首端的无量纲轴向距离X为横轴，以距首端的无量纲径向距离R为纵轴，绘制无气核的情况下，水下高速航行体结构剖面绕流表面压力系数C_p曲线；

其中，P_∞为流场中远方来流的压力；U₀为流场中均匀来流的速度；ρ为流场中流体的密度；P_l是水下高速航行体受到的局部表面压力，根据均匀来流中静止的水下高速航行体表面的压力分布获取；

步骤5：获取无气核的情况下，水下高速航行体结构剖面绕流表面压力系数C_p曲线中压力系数C_p最小位置处的坐标(X₀,R₀)，在该位置附近生成初始半径为R_initial的气泡；

步骤6：令i＝i+1，计算流场的扰动速度势φ_d；

其中，p表示流场控制点；q表示边界积分点；S包括所有流场边界；n为边界法向量，指向流场外为其正向；λ为p点观察流场的立体角；G为格林函数；

步骤7：计算气泡溃灭时产生的射流速度u_i和第i次计算中气泡的半径R_i；

u_i＝▽(U₀X₀)+▽φ_d

其中，Δt_i为第i次计算与第i+1次计算之间的时间步长；S_b为气泡的表面积；

步骤8：将计算得到的气泡溃灭时产生的射流速度u_i量纲化，判断气泡溃灭时产生射流的速度和方向；若射流速度大于20m/s，并且方向指向水下高速航行体结构，则判定产生严重的空蚀效应，水下高速航行体的抗空蚀能力不合格；

步骤9：判断是否达到设定的最大计算次数或者已经产生射流抨击，若未达到设定的最大计算次数且未产生射流抨击，则返回步骤6；否则，执行步骤10；

步骤10：获取计算结果中最大的气泡半径；若最大气泡半径小于2.5，则判定水下高速航行体的抗空化能力不合格；若最大气泡半径小于2.5，则判定水下高速航行体的抗空化能力合格，且最大气泡半径越小，抗空化性能越好；

步骤11：输出水下高速航行体的抗空化空蚀能力评估结果。

本发明的理论研究内容如下：

1、采用基于势流理论的边界元法对水下航行体表面压力进行计算的方法。

2、采用直接边界元法计算水下高速航行体附近微小气核的空化形成特性以及游移空化气泡非球状脉动特性。

3、采用辅助函数法和间接边界元法计算高速航行体和空化气泡周围的速度场和压力场。

4、计算结果可视化方法。

预报水下高速航行体引起的成核空化的三维数值方法：

(1)对高速航行体与气核几何尺度进行无量纲化，取一组特征量对各物理量包括长度、压力、速度、速度势和时间进行无量纲处理，建立具有普遍意义的无量纲化系统和模型。

航行体和气核几何尺寸，使用气核的初始半径R_initial进行无量纲化，假定初始气核为球形且初始时刻处于无法向膨胀速度状态。特征压力取为ΔP＝P_∞-P_c，其中P_∞、P_c分别为远方均匀来流压力和可冷凝气体饱和蒸汽压，本发明未考虑重力的影响，因此忽略P_∞在重力方向的改变。以R_initial(ρ/ΔP)^1/2、R_initial(ΔP/ρ)^1/2以及(ΔP/ρ)^1/2作为时间、速度势和速度的特征量。

(2)采用基于势流理论的边界元法对无气核与航行体耦合作用的简单工况进行计算，得到无气泡扰动下绕流航行体的表面压力分布。

无气核扰动的航行体绕流表面压力可采用直接边界元法，将初始气核置于远方，利用非定常伯努利方程近似获得无气核扰动的航行体表面压力。

(3)假定计算初始条件包括：高速航行体结构与气核的尺度比、气核的初始位置，采用三维边界元法求解边界积分方程，计算气核时域运动特性，记录每个时刻气核表面位置与速度，获得气核在运动过程中的迁移和脉动特性，并采用加权最小二乘法对气核表面位置和速度势进行光顺，保证数值稳定性。

数值计算中，假定气核处于高速航行体附近确定无量纲的距离参数、确定航行体与气核的尺度比，进行程序初始化，采用二阶龙格库塔法对气泡表面位置和速度势进行显式更新，采用加权最小二乘法光顺气泡表面位置和速度势，以保证数值稳定性。

(4)根据目标时刻的气核边界和航行体边界的位置、速度和速度势等物理量，采用间接边界元法，求解流场速度；再利用辅助函数法和非定常伯努利方程求解流场压力，避免速度势的差分运算，提高计算精度和稳定性。

气核与航行体周围的压力场和速度场计算采用以下方法，当已知气核与航行体边界速度势后，采用间接边界元法求解边界分布源密度，然后通过分布源密度求解流场速度，再利用辅助函数法和非定常伯努利方程求解流场压力，避免速度势的差分运算，提高计算精度和稳定性。

(5)将三维的气核与高速航行体结构剖分后使用Matlab中的trimesh命令对三维的气核与结构剖面进行显示；将间接边界元的计算结果按照规则的点阵形式输出，通过Matlab中的contourf和quiver命令对气核周围的压力云图和速度矢量进行可视化处理。

进行计算结果后处理中，首先获得气核与航行体对称面一侧节点位置与单元编号，使用Matlab软件中trimesh命令对三维的气核与结构剖面进行显示；压力场与速度场的可视化包括：首先对流场区域进行离散，例如计算气核周围流场中100×100个离散点的速度和压力，使用Matlab软件中的quiver和contourf命令对气核周围流场的的速度矢量和压力进行后处理显示。

本发明针对水下高速航行体附近局部低压区域气核的成核空化进行数值探究，首先根据气核与航行体结构几何尺寸进行无量纲化，建立可以进行类比分析的三维模型，由于实际结构与空化核的尺度比很大，过大的尺度比将会严重影响数值稳定性，并且尺度效应不可忽视，故本发明中选取较小的尺度比以进行计算。针对背景流场中复杂的航行体边界，采用高精度的三维边界元方法进行数值模拟，通过一些数值技术保证数值稳定性，具体体现为气泡周期和半径与轴对称模型计算结果取得良好的吻合。利用间接边界元法求解流场的速度和压力场，有利于揭示气核、结构以及绕流流动三者相互耦合作用的力学机理。因此，本发明实现了一套完整预报水下高速航行体引起的成核空化的三维数值方法。

1、确定气核与结构的尺度比λ，对气核与水下高速航行体结构表面进行离散。

2、气核的初始半径R_initial、气核初生的位置包括轴向位置X与径向位置R影响着气核在航行体流场中的动态行为，故需要假定气核的初始半径和气核初生的位置，假定气核的初始内压P₀等于远方来流压力P_∞，即气核处于初始平衡状态。确定了气泡的初始条件和边界条件后，对程序进行初始化，连续区域流场的控制方程拉普拉斯方程可转化为离散的边界上的边界积分方程，采用二阶龙格库塔法对气泡节点位置r和节点速度势φ更新，由于数值误差的累积，可能导致气泡表面节点扭曲等数值不稳定现象，因此每隔5个时间步对气泡表面节点位置和节点速度势进行光顺。

3、已知流场边界的速度势后，流场内部速度可采用间接边界元法求解。非定常伯努利方程中含有压力项，但

项往往难于精确求解，由于

项同样满足拉普拉斯方程故引入辅助函数，确定辅助函数满足的边界条件后代入边界积分方程求解，相比于传统的差分运算可更精确的求解速度势偏导数项

联合非定常伯努利方程可获得流域内压力场。

4、三维气泡和高速航行体剖面的显示，使用Matlab软件中的trimesh命令处理，其运行格式为：trimesh(element(:,:),x,y,z)。观测区域速度场和压力场的可视化处理，如在对称剖面上气泡周围2L×2L的区域内，L为某一长度，布置n×n个等间距控制点，计算每个控制点的速度和压力，然后使用Matlab软件中的contourf和quiver命令对气泡周围的压力场和速度场进行后处理显示，运行格式分别为：

contourf(x(:,:),y(:,:),p_min:(p_max-p_min)/100:p_max,'Linestyle','None')

quiver(x,y,u,v,'Color','black')

三维气泡和航行体剖面以及该剖面上气泡周围的压力场和速度场画至同一图片，作为最终的后处理效果图。三维气泡剖面图展示了气泡形态的变化，流场中的速度、压力场一定程度上揭示了空化气泡运动特性和空化空蚀形成的机理。

如附图5所示，数值计算程序开始前，将气泡和高速航行体结构表面离散成有限的计算单元，因此需提供高速航行体结构和气核的几何尺寸。高速航行体结构剖面及其绕流表面压力分布如附图1所示，高速航行体是最小压力系数为-0.4的Schiebe头型，X,R分别定义为头型局部距首端的无量纲轴向距离和无量纲径向距离；压力系数C_p定义为

其中U₀为自由来流速度，ρ为水的密度，P_l是局部表面压力，头型表面压力分布可将气核置于远方，使用直接边界元法和辅助函数法近似求解。

水下航行体高速前进的过程中，水中气核的空间分布是任意的，且其尺度可能具有随意性。在数值计算上需假定气核的初始大小R_initial和初始位置(X,R)方能对问题进行合理简化并求解，本发明将初始气泡置于航行体最低绕流表面压力区域附近。为保证数值求解的稳定性，减小尺度效应，航行体与气核的尺度比不能过小。水下高速航行体单方向平动的过程在理论上与静止的航行体在自由来流中绕流的问题等价，故将问题转化为气核随均匀流流过静止的航行体，可将速度势φ分解为均匀流速度势φ_u和扰动速度势φ_d如下：

φ＝φ_u+φ_d＝U₀x+φ_d (1)

其中φ_d可由边界积分方程计算：

式中p表示流场控制点，q表示边界积分点，S包括所有流场边界，n为边界法向量，指向流场外为其正向，λ为p点观察流场的立体角，G为格林函数。为构成定解问题仍需给定速度势边界条件和速度势初始条件，在航行体表面满足不可穿透条件，初始气泡假定无膨胀：

结构表面 (3)

初始气泡表面 (4)

将(1)式代入非定常伯努利方程可得到气泡表面扰动速度势所满足的动力学边界条件：

其中气泡内部的压力P_b可表示为:

(6)式中P_∞表示z平面无穷远处压力即远方来流压力，ξ为表面张力系数，V₀、V分别为初始气泡体积和当前气泡体积，ξ表面张力系数，κ为气泡表面局部曲率，g表示重力加速度。

当外部压力小于P_c时，气体将发生空化效应，若此时气泡内压刚好达到P_c，达到压力平衡的临界状态。

本发明中设定P₀＝1,ζ＝1.4,ξ＝0.0747.，由上式计算得到临界的气泡半径为2.5。若三维数值计算中气泡半径大于2.5证明气核周围流体压力低于饱和蒸汽压，将产生空化。

气泡表面满足的运动学边界条件为：

利用气泡表面满足的动力学边界条件及运动学边界条件，采用二阶龙格库塔法对气泡节点位置矢量r和节点速度势φ进行显示更新。

第i步时的气泡等效半径R_i可由下式计算：

其中Δt_n为第n步的时间步长，S_b为气泡表面。

在时域模拟气泡运动过程中，由于数值误差累计，可能导致气泡表面节点扭曲等不稳定现象，因此每隔5个时间步可采用最小二乘法(三维问题)或五点光顺法(轴对称问题)对气泡表面位置和速度势进行光顺。如附图2所示为λ＝2，初始气核位置为(-5,0)的轴对称与三维空化气泡半径变化对比曲线图，三维数值方法结果与轴对称模型结果取得了很好的吻合，证明了本文提出的预报水下高速航行体引起的成核空化三维数值方法的有效性和准确性。对于非轴对称工况，仅能用三维模型进行计算，如附图三所示为λ＝50，初始气核位置为(48,60)的三维工况气泡(a)膨胀与(b)收缩数值结果图，各簇曲线为对称剖面上的气泡轮廓。

当已知气泡和航行体表面速度势后，采用间接边界元法求解边界上分布源密度，间接边界元法中，流场中任意一点的速度势可由分布源表示：

式中σ表示流场边界上的分布源密度，当已知气泡和航行体表面速度势后，将控制点布置于气泡与航行体边界上，通过上式即可求解出分布源密度σ。然后将控制点置于流场中，进而通过下式求解流场内部的速度：

然后引入辅助函数：

辅助函数在气泡表面与航行体表面满足的边界条件分别为：

通过(11)式即可获得辅助函数ζ对应的分布源密度，最后将控制点置于流场考察点，再次通过(11)式获得流场中的ζ值，进而通过非定常伯努利方程获得气泡在流场中的诱导压力。该方法不同于传统方法采用差分运算计算流场中的

值，提高了计算精度和稳定性。

三维气泡和高速航行体剖面的显示，使用Matlab软件中的trimesh命令处理，其运行格式为：trimesh(element(:,:),x,y,z)。观测区域速度场和压力场的可视化处理，如在对称剖面上气泡周围2L×2L的区域内，L为某一长度，布置n×n个等间距控制点，计算每个控制点的速度和压力，然后使用Matlab软件中的contourf和quiver命令对气泡周围的压力场和速度场进行后处理显示，运行格式分别为：

contourf(x(:,:),y(:,:),p_min:(p_max-p_min)/100:p_max,'Linestyle','None')

quiver(x,y,u,v,'Color','black')

将三维气泡和航行体剖面以及该剖面上气泡周围的压力场和速度场画至同一图片，作为最终的后处理效果图。最后，将计算结果按时间历程进行对比分析，为研究气核在水下高速航行体运动特性提供参考。如附图4所示为气泡坍塌阶段某典型时刻的(a)压力场与(b)速度场。

以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (1)

1.一种水下高速航行体抗空化空蚀能力评估方法，其特征在于，包括以下步骤：

步骤1：输入待评估的水下高速航行体，设定气核的初始半径R_initial，对水下高速航行体的几何尺度与气核的初始半径R_initial进行无量纲化处理，取一组特征量对压力、速度、速度势、时间和密度进行无量纲处理，建立具有普遍意义的无量纲化系统和模型；

步骤2：建立已知边界的均匀来流流场，将无量纲处理后的水下高速航行体模型置于该流场中；设定计数器i和最大计算次数，初始化i＝0；

步骤3：获取当前时刻无气核的情况下，均匀来流流场中静止的水下高速航行体表面的压力分布；

步骤4：以水下高速航行体结构剖面中首端位置为原点，以距首端的无量纲轴向距离X为横轴，以距首端的无量纲径向距离R为纵轴，绘制无气核的情况下，水下高速航行体结构剖面绕流表面压力系数C_p曲线；

其中，P_∞为流场中远方来流的压力；U₀为流场中均匀来流的速度；ρ为流场中流体的密度；P_l是水下高速航行体受到的局部表面压力，根据均匀来流中静止的水下高速航行体表面的压力分布获取；

步骤5：获取无气核的情况下，水下高速航行体结构剖面绕流表面压力系数C_p曲线中压力系数C_p最小位置处的坐标(X₀,R₀)，在该位置附近生成初始半径为R_initial的气泡；

步骤6：令i＝i+1，计算流场的扰动速度势φ_d；

其中，p表示流场控制点；q表示边界积分点；S包括所有流场边界；n为边界法向量，指向流场外为其正向；λ为p点观察流场的立体角；G为格林函数；

步骤7：计算气泡溃灭时产生的射流速度u_i和第i次计算中气泡的半径R_i；

其中，Δt_i为第i次计算与第i+1次计算之间的时间步长；S_b为气泡的表面积；

步骤8：将计算得到的气泡溃灭时产生的射流速度u_i量纲化，判断气泡溃灭时产生射流的速度和方向；若射流速度大于20m/s，并且方向指向水下高速航行体结构，则判定产生严重的空蚀效应，水下高速航行体的抗空蚀能力不合格；

步骤9：判断是否达到设定的最大计算次数或者已经产生射流抨击，若未达到设定的最大计算次数且未产生射流抨击，则返回步骤6；否则，执行步骤10；

步骤10：获取计算结果中最大的气泡半径；若最大气泡半径小于2.5，则判定水下高速航行体的抗空化能力不合格；若最大气泡半径小于2.5，则判定水下高速航行体的抗空化能力合格，且最大气泡半径越小，抗空化性能越好；

步骤11：输出水下高速航行体的抗空化空蚀能力评估结果。

Images