水下连发超空泡射弹的流动与阻力特性研究

2020-11-04施红辉周东辉贾会霞

空气动力学学报 2020年4期

施红辉,周东辉,周栋,贾会霞

(浙江理工大学机械与自动控制学院,杭州 310018)

0 引言

水下连发射弹的超空泡流动问题,来源于水下枪炮的连续发射、超空泡射弹的饱和攻击等海战背景,具有很强的军事应用价值。超空泡射弹是一种新型水下高速武器[1],相比传统的水下作战武器,它显著地提高了自身的运动速度和射程。超空泡射弹的一个重要用途是作为防御武器装备在水面舰船或者潜艇上[2],通过连续发射高速且密集的射弹形成弹幕保护舰船或者潜艇免受鱼雷的威胁。当超空泡武器串联运行时,这就需要考虑2个甚至多个超空泡的相互作用的问题,由于连发超空泡射弹周围流场相互影响、干扰,引起射弹间流场的变化,与水下单发射弹的超空泡有较大的区别。它既包含单发射弹超空泡的特性,又有自身特殊的复杂性。因此,研究水下连发射弹的超空泡演化规律与阻力特性在学术研究和工程应用上都是必要的,这方面的研究还很少见。

俄罗斯和乌克兰的科学家在射弹空泡特性和水动力学方面进行了大量的基础研究[3-7]。Truscott[8]利用高速摄像技术研究了弹体以小攻角入水时弹体表面诱导产生的超空泡流动态特性,并讨论了长径比与弹体头型对入水时产生超空泡的形态影响。Yamashita[9]等对水下超高速的球形与细长体弹体进行了研究,发现超空泡的存在大大降低了水下航行体的阻力。Weiland 和Vlachos[10]通过实验研究了圆柱形射弹高速入水空泡的产生和发展,并指出射弹发射时驱动气体的泄漏促进了超空泡的快速形成。Lee等[11]利用能量守恒原理对射弹高速入水空泡生成、发展过程及闭合特性进行了研究,建立了高速入水条件下的入水空泡动力学模型。国内学者对水下航行体空化问题也进行了一些研究。易文俊[12]等利用Fluent软件对水下高速航行体自然超空化流动进行了数值模拟,获得了航行体的超空泡形态变化特性。金大桥[13]等基于均质平衡多相流理论,利用Fluent软件对水下射弹自然超空泡减阻特性进行了数值模拟,研究了空化数对水下射弹空泡闭合部位和阻力系数的影响,重点分析了水下射弹结构参数对自然超空泡减阻特性的影响,得到了空化器直径、模型长细比和不同尾部形状对水下射弹超空泡减阻特性影响的规律。郭子涛和张伟等[14-15]对几种柱形弹体水平入水进行了实验和理论研究,建立了关于空泡扩展的理论模型,获得了弹体阻力系数和平头、球头以及几种卵形弹体头型系数之间的关系。魏应杰等[16]对超空泡水下航行体的控制问题进行研究,建立了超空泡航行体非线性动力学模型,设计了非线性切换控制策略。施连会等[17]对超空泡航行体的动力稳定性问题进行了数值研究,计算了超空泡航行体的动力不稳定区域和相应的临界频率曲线。鲁传敬课题组[18-19]基于雷诺平均Navier-Stokes方程和混合多相流模型,对通气空泡与自然空泡共同作用下的三维轴对称的通气空泡流动进行了数值模拟,分析了通气量对水动力学的特性影响,给出适当的通气量利于超空泡稳定的结论。施红辉[20-21]课题组利用自行研制的水平超空泡发射系统对水下单发射弹进行了大量试验,观察到了多种工况下的超空泡形态,研究了水下射弹的水动力学特性以及水深和弹体长径比对超空泡弹体阻力系数及空泡形状影响。

从超空泡问题的提出至今,国内外已经有了许多研究成果,但它们基本上都是针对单个水下物体或射弹的。在水下发射导弹或射弹时,武器能否连续发射直接关系到装备的应用和作战效果,这就需要考虑2个甚至多个超空泡的相互作用的问题,而这方面的研究还很少见。何春涛等[22]对在重力加速下低速两射弹串联倾斜入水问题进行了初步试验,发现当两射弹串联间距小到一定程度时,后面的射弹完全进入到前面射弹形成的空泡内部。施红辉等[23]对连发射弹以100 m/s的速度在水下做水平匀速运动诱导的超空泡流场进行了数值模拟,研究发现连发射弹经过初生、发展和相互作用后,能够形成一个包裹两射弹的超空泡,超空泡的尾部产生高速回射流。

针对水下连发射弹自然减速运动过程中的超空泡演化规律及运动特性,本文以Flunet流体仿真软件为平台,利用动网格技术分别建立了水下两连发、三连发射弹的二维数值计算模型,通过数值计算分别获得了两连发射弹和三连发射弹在水下自由运动过程中的超空泡演化规律、连发射弹头部前方的压力分布、连发射弹的位移变化曲线、连发射弹的速度衰减曲线、连发射弹的阻力特性曲线,本文的研究结果率先揭示了水下两连发和三连发射弹诱导的超空泡流动的典型特征及相关力学参数。

1 控制方程及数值方法

1.1 控制方程

数值计算采用流体体积函数(VOF)多相流模型描述水和水蒸汽构成的多相流动系统。VOF多相流模型将水、水蒸汽两相当做单一介质的混合流动系统,各相共享一套动量方程,通过计算得到单元内各相的体积分数确定流动系统中各相的分布。

混合相的连续性方程和动量方程分别为:

式中:t为时间;ui、uj分别为笛卡尔坐标系中的速度分量;ρm为混合相密度;p 为流场的压力;μm为混合相动力黏度;SM为附加的源相。

ρm,μm的表达式分别为:

式中:αv为水蒸汽相的体积分数,ρv、ρl 分别为水蒸汽和水的密度;μv、μl分别为水蒸汽和水的动力黏度。

数值计算采用RNG k-ε 湍流模型[24],由于空化现象是流线强烈弯曲导致的,该模型能够更好地处理高弯曲流线和高应变率的流动。湍动能k 和耗散率ε 的控制方程分别为:

式中:μeff＝μm+μt,μt 为湍流黏度;k 为湍动能,ε 为耗散率,αk、αε分别为k 和ε 的负向效应的普朗特数;Gk为平均速度梯度产生的湍动能;C1ε、C2ε为湍流动能耗散率的经验常数。

1.2 空化模型

空化模型选择Schnerr-Sauer空化模型[25],该空化模型具有较高收敛速度和计算稳定性。水蒸汽相方程的一般形式为:

式中:vv为水蒸汽相的速度矢量;Re、Rc分别为蒸发速率和冷凝速率;ℜB为气核的半径;Pv为水的饱和蒸汽压力。

1.3 数值方法和边界条件

连发射弹的计算模型为轴向串联的两发射弹和轴向串联的三发射弹,定义相邻两射弹的间距为H,射弹的物理模型是直径D＝6 mm,长度L＝48 mm的圆柱体,如图1所示,来源于参考文献[21]的试验模型,材质为5005铝镁合金,质量为3.56 g。由于数值模拟中射弹设定的速度很高,计算求得的弗劳德数Fr≫30,根据重力效应的判断准则,可以不考虑重力的影响。由于射弹为轴对称回转体,数值计算采用二维轴对称模型。整个网格区域宽度是200 mm,长度是1000 mm,足够避免边界效应。数值模拟中采用射弹运动、流体静止的方法研究连发射弹的自然超空泡,这种方法更接近真实的环境。通过用户自定义函数(UDF)编写程序嵌入到Flunet中模拟具有初动能的连发射弹自由运动,在运动过程中网格的更新和生成采用动态层法实现。分别对水下两连发和三连发射弹进行了数值模拟,连发射弹彼此的间距H ＝10D＝60 mm,这两种情况的计算模型如图2所示。图2(a)表示两连发射弹的计算域及边界条件设置,图2(b)表示三连发射弹的计算域及边界条件设置,这两种计算模型的计算域边界条件划分相同,右侧边界条件为压力出口边界条件,左侧边界条件为无滑移壁面条件,计算域上侧设定为无滑移壁面条件,射弹设置为无滑移壁面条件,射弹从计算域右侧运动到左侧,射弹的初始发射速度都设定为相同的速度v0＝200 m/s,之后射弹做自由运动。整个计算域采用四边形结构网格划分,射弹附近的网格进行加密处理。

图1 射弹的物理模型Fig.1 Physical model of projectile

图2 连发射弹的计算域和边界条件设置Fig.2 Computational region and boundary conditions of successively fired projectiles

对于上述计算模型压力与速度之间的耦合求解采用PISO 算法,压力场的空间离散采用PRESTO!格式,各项体积率离散采用CICSAM 格式,密度和动量采用二阶迎风离散格式,空间离散采用二阶迎风格式,时间离散采用一阶隐格式。

2 数值方法验证

在研究水下连发射弹超空泡流动前,先进行数值方法的验证。由于水下连发超空泡射弹试验难度很大,并且缺乏相关的文献数据,所以选取水下单发超空泡射弹试验进行数值模拟验证。采用上文所述的数值模拟方法对参考文献[21]中的工况8进行了数值计算,工况8采用的射弹和图1所示的射弹模型相同,质量为3.56 g,试验中射弹的初速度大小为83.1 m/s,水深为90 mm。数值计算的整个网格区域宽度是200 mm,长度是1000 mm,右侧边界条件为压力出口边界条件,左侧边界条件为无滑移壁面条件,网格为四边形网格,射弹周围的网格进行加密处理,加密区的网格尺寸为0.2 mm,时间步长取4×10－7s。

在水下深度一定的情况下,射弹高速无推力自由运动过程中,忽略重力和浮力的影响,射弹形成的超空泡内外压力差Δp 相同,空化数仅由射弹速度vp决定,空化数σ 的表达式[26]如下:

式中:h 为水的深度,p0为一个标准大气压,pc为水的饱和蒸汽压,这里取值为3540 Pa。

超空泡状态下圆柱体的阻力系数可以近似表达为[23]:

式中:Cd0是空化数σ＝0的阻力系数[26](圆盘空化器的Cd0≈0.82)。

对于水下单发射弹无推力自由运动,由牛顿第二定律可得:

式中:m 是射弹质量,S 是射弹位移,A0是射弹的截面积。

定义速度衰减系数kv＝ρlA0Cd0/2m ,通过积分方程(12),可得:

关于单发超空泡射弹的位移随时间变化,图3显示了数值模拟结果、实验结果、公式(13)的理论计算结果三者之间的对比,纵坐标表示的是无量纲位移(S/D),从图中可以看出数值模拟的射弹位移变化趋势与实验和理论计算结果基本一致,在数值上吻合良好。图4比较了数值模拟和实验在t＝2 ms时获得的射弹水下空泡形状,从图中可以看出两者的空泡形状基本一致。经过上述的对比,有效地验证了数值模拟方法的准确性。

图4 数值模拟和实验在t＝2 ms时的超空泡形状Fig.4 Comparisons of supercavity configurations between numerical simulation and experiments at t＝2 ms

以水下两连发射弹为例进行网格无关性验证,建立了3 中不同密度的网格,网格数分别为45.5 万(case1,网格最小尺寸为0.4 mm)、53.7万(case2,网格最小尺寸为0.2 mm)、65.8万(case3,网格最小尺寸为0.1 mm),取时间步长为4×10－7s,分别采用上述三种网格密度对发射速度为200 m/s的两连发射弹入水进行了数值模拟,得到了前发射弹速度随时间的变化曲线,如图5所示。从图中可知,随着网格密度的增加,case2的网格密度和case3的网格密度计算所得的前发射弹运动速度之间的差异已经小到可以忽略不计,综合考虑计算精度和计算成本,选取了网格密度为case2的网格用于本文计算。

在不影响计算正常的前提下,选取了三种时间步长值(d t1＝4×10－7s,d t2＝2×10－7s,d t3＝1×10－7s),针对不同的时间步长对两连发射弹入水的计算结果进行了分析验证,获取了前发射弹速度随时间的变化曲线,如图6所示。通过对比可知,三种时间步长下前发射弹的速度变化规律一致,这说明此范围内的时间步长对于本文的数值模拟是适用的,考虑时间成本的限制,选择了时间步长d t1＝4×10－7s进行数值模拟。

图5 不同网格密度的计算结果Fig.5 Calculated results of different mesh densities

图6 不同时间步长的计算结果Fig.6 Calculated results of different time steps

3 结果与讨论

水下连发射弹诱导的自然超空泡的发展是一个非定常过程。水下两连发和三连发射弹的超空泡演化过程如图7所示,图7中红色代表水相,蓝色代表水蒸汽相。图7(a)表示的是水下两连发射弹超空泡演化过程,对其进行分析,可以将两连发射弹超空泡的发展过程分成3个阶段:第1阶段,射弹1和射弹2各自形成独立的超空泡,射弹被包裹在空泡内部,只有射弹头部沾湿,射弹受到的阻力主要为压差阻力。第2阶段,两射弹运动过程中,在射弹1超空泡流场的影响下,射弹2的头部逐渐靠近并最终进入射弹1超空泡的尾部,两发射弹的超空泡发生部分融合,随后第二个超空泡的大部分被分离(图中对应t＝1.0～1.2 ms),分离后的空泡在周围水压的作用下慢慢变小并最终发生溃灭,射弹2进入射弹1的超空泡内部运动后并逐渐追上射弹1。t＝1.5 ms时,射弹2已经完全在射弹1超空泡的内部运动。第3阶段,射弹2与射弹1发生追尾碰撞,碰撞后两者共同运动,文献[22]在开展的两个平头圆柱体串联低速入水试验中有观察此种现象,这也证明了数值模拟结果的准确性。t＝2.6 ms时,射弹2追上射弹1并发生追尾,在两射弹发生追尾后,超空泡壁面上出现一个波动现象,射弹1头部流动分离点附近的空泡壁面先出现波谷(图中t＝2.6 ms时,①处的放大图中椭圆标记处),之后发生追尾处的空泡壁面出现波峰(图中t＝2.7 ms时,②处的放大图中椭圆标记处),推测此现象产生的原因是由两射弹发生追尾碰撞产生的瞬态扰动引起的;该波动现象一直存在,随着时间的推移,存在波动的空泡壁面出现颈缩,空泡从此位置逐渐开始分裂。图7(b)给出了水下三连发射弹超空泡演化过程,其发展过程可以分为4个阶段:第1阶段,三发射弹各自形成独立的超空泡。第2阶段,射弹3进入射弹2的超空泡,射弹2进入射弹1形成的超空泡。第3阶段,射弹3从射弹2脱离的超空泡内部穿出,穿出时其头部与水接触再次生成空泡,生成的空泡与射弹2脱离的空泡相互融合;第4阶段,射弹2先追赶上射弹1并发生追尾碰撞,碰撞后两者共同运动,随后射弹3追赶上射弹2并发生追尾碰撞,碰撞后三者共同运动。t＝2.9 ms时,射弹2和射弹1发生追尾碰撞,射弹3 部分进入射弹1 的超空泡内。t＝4.9 ms时,射弹3与射弹2发生追尾碰撞。射弹发生碰撞后同样在超空泡壁面上出现扰动现象,与两连发射弹的情形一样,这里不再详细叙述。随着时间的推移,三连发射弹存在扰动的空泡壁面出现颈缩,并且空泡从此位置开始分裂。

图7 连发射弹水下运动诱导的超空泡演化过程Fig.7 Evolution of supercavities induced by successively fired projectiles

图8 表示的是两连发射弹在时间t＝0.7 ms时的压力云图以及射弹头部中心点前4D 处的压力分布,图8(b)中的横坐标为无量纲位移S1/D,其中S1表示点到射弹头部中心点的水平位移,纵坐标表示无量纲压力p/p0。从图中可知,高压区出现在前、后发射弹的头部表面,前发射弹头部中心点处的压力达到85倍的大气压,后发射弹头部中心点处的压力为57倍的大气压,在射弹头部前方流域压力值迅速降低,距离平头顶点越远,压力越小;低压区出现在两发射弹的空泡内,由于空化效应其压力值和饱和蒸汽压力保持一致(约为3540 Pa)。此时刻前发射弹超空泡的尾部与后发射弹的头部相距的距离很小,从而引起两发射弹周围的流场互相干扰,具体表现为前发射弹超空泡尾部空泡在后发射弹头部附近的高压区作用下部分溃灭,后发射弹头部中心点处的压力及头部前方的高压区范围小于前发射弹的原因是前发射弹超空泡尾部存在回射流[19],导致后发射弹相对前方流体的速度减小,而且前发射弹形成的空泡射流为汽水混合物,密度远小于水的密度。

图9给出了三连发射弹在时间t＝0.7 ms时的压力云图以及射弹头部中心点前4D 处的压力分布。从图中可以看出射弹1头部驻点处的压力最大,其大小约为85倍的大气压,而射弹3和射弹2头部前方的压力值比较接近,这说明后发射弹都会受到前一发射弹所形成的超空泡流场的影响。

图8 t＝0.7 ms,两连发射弹的压力分布图Fig.8 Pressure distribution of two successively fired projectiles at t＝0.7 ms

图9 t＝0.7 ms,三连发射弹的压力分布图Fig.9 Pressure distribution of three successively fired projectiles at t＝0.7 ms

图10 给出了连发射弹速度随时间变化的曲线。图10(a)表示的是两连发射弹速度变化,可以看出射弹1和射弹2在运动初期的速度衰减趋势相同,并且速度衰减很快,这是由于两射弹的速度开始很高,受到水的阻力大;在0.7 ms＜t＜2.6 ms,射弹2的速度衰减逐渐变慢到速度保持不变,射弹1的速度继续衰减,这是因为射弹2逐渐接近射弹1尾部空泡时,受到前发射弹尾部空泡低压区的影响,受到的压差阻力减小,当射弹2进入到射弹1的空泡中时,射弹2周围的介质变成了水蒸气相,密度远小于水,导致射弹2受到很小的阻力,这里我们不妨标记为0,而射弹1继续受到水的阻力作用;在时间t＝2.6 ms,两射弹的速度发生跃变,射弹1的速度突然变大,射弹2的速度突然变小,这是因为射弹1和射弹2发生了追尾碰撞。碰撞后,射弹1、射弹2的速度衰减相同。图10(b)表示的是三连发射弹速度随时间的变化。在射弹运动初期,射弹1、射弹2、射弹3 速度的衰减相同。在1 ms＜t＜2.1 ms,射弹3开始进入射弹2的超空泡中,所以此时间段射弹3的速度基本保持不变;在2.1 ms≤t＜2.5 ms,射弹3从前一发射弹脱离的超空泡中穿出,受到水的阻力其速度又开始减小;2.5 ms≤t＜4.9 ms时,射弹3进入了射弹1的超空泡中导致其所受阻力几乎为0,从而速度保持不变。射弹2 在1.2 ms≤t＜2.9 ms进入了射弹1的超空泡中引起它的速度保持不变,t＝2.9 ms时,射弹1和射弹2发生追尾碰撞引起两射弹的速度发生跃变;2.9 ms≤t＜4.9 ms,射弹1和射弹2速度变化相同。t＝4.9 ms时,射弹3与射弹2发生碰撞,射弹1、射弹2、射弹3三者的速度发生跃变;t＞4.9 ms,射弹1、射弹2、射弹3速度衰减相同。

图10 连发射弹速度衰减曲线Fig.10 Velocity attenuation of successively fired projectiles

图11 显示的是连发射弹位移随时间变化的曲线,其中纵坐标为无量纲位移S/D,坐标原点为初始时刻射弹的重心位置,位移的正方向为射弹运动的速度方向。图11(a)表示的是两连发射弹的无量纲位移变化,从图中可以看出随着时间增加,射弹2的位移相比射弹1的位移越来越大,证明两射弹之间的距离在逐渐缩短,在t＝2.6 ms,射弹2的位移恰好比射弹1的位移大10D,说明射弹2在此时追上射弹1,之后两射弹的位移差始终为10D,表明追尾后射弹1和射弹2共同运动。图11(b)给出了三连发射弹无量纲位移的变化,射弹2与射弹1的位移差随时间逐渐增加,射弹3与射弹2的位移差随时间逐渐增加,在t＝2.9 ms时,射弹2 的位移恰好比射弹1 的大10D,表明此时射弹2追上射弹1,t＝4.9 ms时,射弹3的位移恰好比射弹2的大10D ,表明此时射弹3追上射弹2。

图11 连发射弹位移随时间变化曲线Fig.11 Displacement of successively fired projectiles with time

射弹在水中运动时,会受到摩擦阻力和压差阻力的影响,这两部分阻力的和为射弹受到的总的阻力。对该阻力进行无量纲化,则得到无量纲阻力系数。无量纲阻力系数[23]的定义如下:

式中:Fd为射弹阻力。

根据式(14),可得到如图12所示的连发射弹的阻力系数变化曲线图。图12(a)表示的是两连发射弹的阻力系数随时间(射弹1和射弹2发生碰撞前)的变化,在运动的起始阶段,射弹1和射弹2刚与水接触,射弹要带动周围的流体质点一起运动,受到水的阻力很大,所以此阶段射弹1、射弹2的阻力系数较大;随着射弹1与射弹2的超空泡尾部逐渐接近,射弹2受到射弹1超空泡流场的影响,从而引起射弹2头部附近的压力降低,导致射弹2受到的压差阻力减小,对射弹2而言,射弹1相当于起到导流的作用,因此射弹2的阻力系数逐渐降低;随着射弹2进入到前一发射弹的空泡中运动,此时射弹2的阻力系数降低到几乎为0,而射弹1继续受到水的阻力,它阻力系数在某一定值附近震荡变化。

图12(b)表示的是三连发射弹的阻力系数随时间(射弹1和射弹2碰撞前)的变化曲线。射弹1的阻力系数在起始时刻较大,之后逐渐减小到某一定值附近振荡变化;射弹2的阻力系数随时间增加逐渐减小到几乎为零,这是由于射弹2在逐渐接近射弹1的过程中,受到射弹1的超空泡流场的影响所致;射弹3的阻力系数随时间增加先减小到趋近于0,这是因为当射弹3进入射弹2形成的超空泡内部时所受阻力很小,在时间t＝2.1 ms时射弹3从射弹2的超空泡内部穿出,其头部与水突然接触受到一个巨大的瞬时作用力,导致射弹3的阻力系数突然增大到一个很高的值,之后射弹3进入射弹1形成的超空泡内,射弹3的阻力系数又减小到接近于0。