吴霆锋, 杜永成, 杨立, 王保霖, 张泽华
(海军工程大学 动力工程学院,湖北 武汉 430033)
潜艇航行过程中排出大量温度较高的冷却水,冷却水浮升扩散,可能使海面上尾流区域出现明显的冷热特征,将这种尾流称之为热尾流。理论上,基于红外技术对潜艇尾流水面冷热特征进行分析,可以获得潜艇的航向、航速、艇型等至关重要的作战信息。同时,海洋环境中存在盐度环境,即海水盐度可能随深度不同而变化的现象,盐度环境对于冷却水的浮升扩散可能存在一定影响,进而影响海面尾流区域的冷热特征。因此,开展盐度环境下的潜艇热尾流研究有重要意义。
杨立等[1]通过实验探究了稳定温度密度分层环境对于水下航行器热尾流浮升的抑制作用;刘兴泉[2]整理得到了东海P-N断面夏季温盐度的分布特征;Bonnier等[3]对线性稳定盐层中球体尾迹的扩散过程进行了实验研究;Fritts等[4]基于DNS方法模拟了潜艇尾迹演化过程;Meunier等[5]基于流体的自保性,建立了线性分层流体中拖曳尾迹的一般模型;戴天奇等[6]基于Fluent动网格法对潜艇水下航行的二维计算模型进行了仿真,结果表明,动网格方法具有较高精度,能反映尾流浮升过程中的瞬时变化;王平等[7-8]基于重叠网格技术进行了均匀水体中的潜艇热尾流仿真模拟,并与传统来流法进行了对比,发现了重叠网格技术具有更高精度。
综上所述,国内外对于潜艇尾流在海洋环境下的发展过程已经有了一些理论和实验研究,但大多采用传统来流法,未采用更为先进与精确的重叠网格方法;且多数未考虑海洋盐度对尾流的影响研究,缺少系统性研究海洋盐度对于潜艇热尾流的影响。因此,本文将基于重叠网格技术,开展不同盐度环境下的潜艇尾流仿真,采用有限体积法进行模拟计算,探究海洋盐度对于潜艇尾流浮升扩散与水面冷热特征的影响。
对海洋盐度环境中的水下航行潜艇尾流浮升扩散过程进行数值模拟,需满足如下控制方程。
连续性方程:
(1)
动量方程:
(2)
能量方程:
(3)
组分传输方程:
(4)
式中:T为温度,K;ρ为流体密度,kg/m3;ui、uj(i、j=1,2,3)为x、z、y方向的速度分量,m/s;gi为重力加速度分量,m/s2;μ为粘性系数,Pa·s;λ为导热系数,W/(m·K);cp为定压比热容,J/(kg·K);cs为组分的体积浓度;ρcs为组分的质量浓度;Ds为组分的扩散系数,m2/s。
采用物性多项式构建海水密度与温度的关系[9]:
ρ(T)=A+BT+CT2+DT3+ET4+FT5
(5)
式中各项系数如表1所示。
表1 海水物性多项式参数[9]Table 1 Polynomial parameters of seawater physical properties
以缩比潜艇模型为研究对象,按照1∶100的比例建立计算模型。模型中计算域的长×宽×高为10.32 m×1.0 m×1.126 2 m,分别沿x、z、y方向,艇长为86 cm,艇体直径为10 cm,两侧冷却水排放口直径为10 mm,以初始时潜艇中心为模型原点。潜艇航速为0.1 m/s,冷却水排放速度为0.05 m/s,温度为350 K。冷却水排放口设置为速度入口,艇体周围与背景水域交界面设置为重叠界面,顶部水面出口为压力出口,四周壁面为对称壁面,底部设置为静止壁面。
图1为基于重叠动网格的潜艇尾流三维计算模型示意图。
图1 基于重叠动网格的潜艇尾流三维计算模型示意Fig.1 Schematic diagram of the submarine 3-D calculation model based on overset grid
在Fluent Msehing平台进行模型网格划分。潜艇及近艇区域的重叠部分采用Poly-Hexcore(多面体-六面体混合网格)进行网格划分,背景海水区域采用Cartesian(笛卡尔网格)方法进行结构化网格划分。为进行网格无关性验证,将模型网格数目分别划分为101万、202万、301万、407万,令潜艇进行t=9 s的航行仿真,选取距离潜艇后方0.2 m处截面,截面垂直于潜艇的运动方向。测得截面上热尾流中心位置处的温度,结果如图2所示,可以看出,当网格数目达到202万时,截面上温度不再随网格数的增加而变化,故选定总网格数目为202万。
图2 网格无关性验证Fig.2 Mesh amount independence verification
湍流模型选用Realizablek-ε模型。压力-速度耦合模式选用Coupled算法,压力项选用Body Force Weighted,动量方程、能量方程、湍动能和耗散率均选用二阶迎风格式。采用瞬态计算,时间步长设置为0.01 s,计算时长为60 s。
海洋盐度指的是海水中全部溶解固体与海水重量之比,通常以每千克海水中所含的克数表示。在海水溶解物中,NaCl占据主要部分。因此,可将盐度近似表示为海水中盐类物质的质量分数。据资料显示,世界大洋的平均盐度为35‰。
一方面,海洋盐度受纬度影响,因海域所处的纬度位置不同而具有不同的盐度;另一方面,盐度还主要受到蒸发量和降水量之差的影响。水分蒸发使得海水浓缩,从而盐度上升,降水使海水稀释,盐度下降。在本文中,计算域尺寸为10.32 m×1.0 m×1.126 2 m,纬度变化较小,因此主要考虑蒸发量和降水量之差的影响。
当海水蒸发量大于降水量时,海水上层盐度上升,下层海水盐度维持稳定,故呈现出上层盐度高、下层盐度低的正梯度盐度分层分布,此类正梯度盐度分层属于不稳定分层,海水会受重力作用产生上下对流,分层存在时间短暂,故不予考虑。当海水蒸发量小于降水量时,海水上层盐度下降,下层海水盐度变化缓慢,基本维持原有盐度,因此呈现出上层盐度低、下层盐度高的负梯度盐度分层分布,目前海洋环境中以该类分布为主,这也是本文主要探讨的盐度分层环境。
研究盐度分层环境时,简化为海水盐度沿深度方向线性变化,当盐度分层的上下边界密度均给定时,即可计算得出表征分层流场特征的浮频率数N,其计算式为:
(6)
式中:ρ0为上边界流体密度,kg/m3;Δρ为上下边界流体密度差,kg/m3;H为流场深度,m。浮频率数作为分层流场的特征数,当其较大时,代表流场密度分层跨度较大,密度分层梯度大;反之,当浮频率数较小时,代表流场的密度分层梯度较小。
因NaCl溶质为固体,其物性与液态水存在较大差距,模拟其混合溶解过程存在较大难度,且误差较大。因此,考虑将NaCl溶液作为替代性介质,引入体积数20%NaCl溶液,其主要物性参数如下[10]。
20%NaCl溶液的密度多项式为:
ρNaCl(20%)=-0.001 3T2+0.28T+1 176.37
20%NaCl溶液的粘度多项式为:
μNaCl(20%)=-4.80×10-9T3+4.88×10-6T2-
1.66×10-3T+0.19
20%NaCl溶液的导热系数多项式为:
kNaCl(20%)=-5.75×10-6T2+
4.81×10-3T-0.33
20%NaCl溶液的比热容多项式为:
cpNaCl(20%)=0.94T+3148.68
在Fluent中引入体积分数20%NaCl溶液后,通过初始化Ptach功能设置其分布,模拟形成海洋盐度环境,具体设置情况见表2。
表2 盐度分层参数设置Table 2 Parameter of temperature stratification
本节将对均匀盐度与盐度分层2类海水环境数值仿真的结果进行分析。一方面,海面温度特征是潜艇热尾流探测最直接的依据,明显的海面温度特征容易被红外探测捕捉到,进而暴露潜艇踪迹,因此,潜艇航行时的海面温度场是结果分析的重点。另一方面,对潜艇冷却水浮升扩散过程的研究,能更深刻地认识潜艇热尾流的形成机理。
因此,可分别通过对海水温度场与冷却水浮升扩散过程进行分析,研究盐度对于潜艇尾流的影响。
2.1.1 海面温度特征
选取上界面温度场进行分析,图3为t=60 s时不同盐度的均匀海水环境下的海洋表面温度云图。观察图可以看到,不同盐度下均匀海水环境的温度分布基本一致,合理推断均匀海水环境下,环境盐度的改变对潜艇尾流形成的海洋表面温度场影响非常微弱。
图3 均匀盐度下的海洋表面温度云图Fig.3 Temperature cloud of sea surface in uniform sea water
2.1.2 冷却水浮升扩散
在研究潜艇冷却水的浮升扩散过程中,冷却水中心的浮升高度与温度变化是重要指标。冷却水中心浮升高度即冷却水中心的y坐标,可以直观反映冷却水的浮升轨迹;冷却水中心的温度变化则能清楚反映冷却水与海水的热交换程度。
图4为潜艇冷却水中心温度随航行距离变化图,图5为冷却水中心浮升高度(即冷却水中心的y坐标)随航行距离变化的规律。观察图4,可以发现,3种不同均匀盐度环境下冷却水温度变化与浮升过程基本一致,均匀海水下整体盐度的改变对于潜艇冷却水的浮升扩散影响微弱,与上述温度场分布结果吻合,与上述结论相互印证。需要说明的是,在图5中出现冷却水中心短暂下降的现象,该时段为冷却水浮升至潜艇艉翼附近,向上浮升受到艉翼阻挡后反弹,故而短暂下沉。
图4 均匀盐度环境下冷却水中心温度Fig.4 Central temperature of cooling water in uniform seawater
图5 均匀盐度环境下冷却水中心浮升高度Fig.5 Buoyancy height of cooling water center in uniform seawater
2.2.1 海面温度特征
图6所示为t=60 s时的负盐度梯度分层环境下的海面温度云图。观察图6可以看到,在盐度分层环境下,当盐度梯度增大时,浮频率随之增大,海面尾迹区域明显缩小,温度降低。因此可以得出结论,盐度分层环境下海面温度特征随盐度梯度的增大而减弱。
图6 盐度分层环境下的海面温度云图Fig.6 Temperature cloud of sea surface in salinity stratification environment
为了更清晰地研究海面温度变化,选取海面上直线AB与直线CD,进行两直线上的温度分布分析,直线AB、CD分布如图6所示。图7、图8分别为直线AB与CD上的温度分布。观察发现,当盐度分层梯度增大,即浮频率增大时,海面温度分布趋势保持不变,尾流区域整体温度降低。因此可以得出结论:盐度分层梯度的增大对潜艇热尾流形成的海面温度场具有一定的抑制作用。
图7 直线AB温度分布Fig.7 Temperature distribution on line AB
图8 直线CD温度分布Fig.8 Temperature distribution on line CD
图9所示为海面最高温度随浮频率数的变化规律,可以发现,浮频率数N=0的3种均匀盐度环境的海面最大温度基本相同,最大温度均处于291.024 K左右,而当盐度存在分层,且分层梯度产生变化时,随着浮频率数N逐渐增大,海面最大温度先略微上升,达到291.026 K后降低,降低速度随N的增大而减小,整体变化趋势接近正弦函数。
图9 海面最高温度随浮频率数N变化规律Fig.9 Variation of maximum sea surface temperature with buoyancy frequency N
2.2.2 冷却水浮升扩散
图10所示为t=60 s时不同盐度分布下的中性面温度云图。分析图10可以发现,在盐度分层环境下,潜艇冷却水浮升受到明显抑制。当浮频率数增大时,冷却水被释放后向后方扩散的趋势加剧,向上浮升的趋势受阻。
图10 盐度分层下的中性面温度云图Fig.10 Temperature cloud of neutral surface in salinity stratification environment
图11所示为潜艇冷却水中心浮升高度,可以发现,随着浮频率数的增大,潜艇冷却水浮升受到明显抑制,最大浮升高度有所降低,甚至出现到达浮升最高点后下沉的现象。图12为盐度分层环境下潜艇冷却水中心温度随航行距离的变化规律。可以看出,当浮频率较小时,冷却水中心温度较低,而冷却水排放之初温度均为350 K,说明浮频率数较低时,冷却水与周围海水的热交换更加剧烈,而当浮频率数增大时,热交换减弱,从而导致尾流浮升受到抑制,进而使得海洋表面的尾流区域温度整体降低。
图11 盐度分层下冷却水中心浮升高度Fig.11 Central buoyancy height of cooling water under salinity stratification
图12 盐度分层下冷却水中心温度Fig.12 Central temperature of cooling water under salinity stratification
2.2.3 旋涡结构
潜艇热尾流问题本质上接近于横流环境下的热射流排放问题,研究潜艇尾流所引起的涡结构,有助于更深入认识潜艇热尾流的发展演变规律。图13所示为N=0.076 87时的盐度分层下的潜艇后方yz平面速度矢量图,位置为x=-5.4 m。观察发现,艇体上下均存在一对反向旋涡。潜艇两侧排水口向外排放冷却水时,冷却水受排放惯性作用向下运动,便产生了艇体下方的旋涡对;继而高温冷却水受到浮升力作用向上浮升,在艇体上方形成了另一对旋涡。
图13 N=0.076 87时潜艇后方yz平面速度矢量图Fig.13 Velocity vector diagram of yz plane behind submarine when N=0.076 87
图14所示为上下旋涡对的旋涡半径随浮频率数的变化规律,可以发现,旋涡半径随浮频率增大而减小,减小速度随浮频率增大而逐渐降低,且艇体下方旋涡半径的减小速度比上方旋涡更快。因此,盐度分层梯度的增大对于潜艇尾流所引起的旋涡具有抑制作用,促使其影响面积收缩,抑制作用随着盐度分层梯度的增大逐渐减弱。同时,盐度分层对于艇体下方的旋涡抑制作用更为明显。
图14 上下旋涡对半径随浮频率数N的变化规律Fig.14 Variation of vortex radius of upper and lower vortex pairs with buoyancy frequency N
1)均匀海水环境下,环境盐度的改变对冷却水的浮升扩散影响非常微弱。
2)盐度分层环境对潜艇冷却水浮升具有明显抑制作用,海面温度特征随盐度梯度的增大而减弱,尾流区域最高温度与浮频率基本呈正弦函数关系。
3)盐度分层对于潜艇航行产生的反向旋涡具有抑制作用,抑制作用随着盐度分层梯度的增大逐渐减弱。