0 引　言

潜艇，特别是续航能力强的核潜艇，在海洋中航行时，由于需要冷却动力装置，因此会向海洋中排放大量的冷却水而形成热尾流^[1]。在浮升力作用下热尾流漂浮至海面，引起海面温度分布产生变化。通过分析海面温度特征，理论上可以获得潜艇航向、位置等重要信息，从而实现对潜艇的跟踪与打击。因此，开展潜艇热尾流浮升规律的研究具有重要意义。

早在1972年Merrit等^[2]便开始了对潜艇热尾流温度分布特性的研究，他们利用热敏电阻对模拟的潜艇无动量尾流进行探测。1988年，俄罗斯科研人员，研究指出对于潜艇而言，热尾流红外探测是一种有效的非声探测方法^[3]。之后国内外学者又进行了大量的实验研究，并取得了一系列的重要研究成果^{[4 – 7]}。但热尾流浮升扩散是一个极为复杂的过程，实验研究成本高、周期长，且易受测量方法及仪器精度的限制。随着计算机技术的发展，数值方法以其低成本、高效率的优势，而越来越成为一种重要的科研手段。为方便研究，早期研究人员大都忽略海面空气流动的影响，将海面简化为壁面或压力边界。但真实环境中海面是大气和海洋的交界面，海面上空气流动、交界面内气液两相能量、动量交换都会对研究结果产生重要影响。吴猛猛等^[8]和张昊春等^[1]采用VOF法在二维空间内对潜艇热尾流浮升扩散过程进行仿真研究，指出在水下热尾流浮升主要受浮力和海水垂直温度梯度的影响。到达水面后热尾流将与气液交界面发生掺混，逐渐扩散和衰减。研究均采用相对运动法模拟潜艇航行，即将潜艇模型视为静止，根据相对运动原理，让海水以与航速相同的速度向模型流动。这与潜艇航行的实际情况仍有一定差距。戴天奇等^[9]基于动网格技术，在二维平面内对多种工况下潜艇冷却水排放过程进行仿真。研究表明动网格技术对于潜艇航行过程中热尾流浮升过程有较高的计算精度。Shi Zhi-guang等^[10]采用动网格技术，对不考虑海面空气流动条件下潜艇热尾流海面温度分布特征及其红外识别过程进行仿真分析，并获得了热尾流水面温度分布及其红外成像特征。

综上所述，国内外学者虽然对热尾流浮升规律已经进行了广泛的研究，但在三维空间内考虑海面空气流动影响的热尾流水面温度分布特征方面的研究尚不多见。本文基于动网格技术，分别采用壁面法与VOF法对三维海域内潜艇热尾流浮升过程进行仿真研究。获得了潜艇热尾流的浮升、扩散规律及其水面温度分布与红外识别特征，分析了水面空气流动及大气-海水交界面对热尾流发展过程的影响，并通过对不同冷却水排放方案下热尾流水面温度分布特征的变化规律的分析，提出了抑制热尾流浮升的方法。本文研究结果可以为潜艇红外探测及红外隐身方面研究提供理论基础。

1 数理模型 1.1 计算模型

图1为本文研究模型示意图，模型长1 200 m，水面宽200 m，水深60 m（在VOF法中增加10 m空气层，模型总高度为70 m）。潜艇在水下40 m深度处沿海水流动方向航行。

本文以 “海狼”级核潜艇为研究对象，表1为其基本参数，研究重点在于潜艇热尾流在水下浮升过程及其在水面温度分布特征。为方便计算，文中忽略了潜艇瞭望塔及方向舵等结构，将其简化为长110 m，最大直径10 m，排水口直径2 m的长柱体结构，如图1所示。

1.2 网格技术

图2为计算模型网格局部示意图，为提高计算精度，对潜艇周围及气液交界面处网格进行了局部加密处理。

网格划分质量及网格数量对求解结果的精度有重要的影响。为保证本文模拟结果的精准性，文中网格质量均大于0.4，角度质量均大于36°，最大纵横比均小于30。为验证网格无关性，本文分别划分了278万、516万、900万、1 022万和1 796万5种数量的网格，分别在距离潜艇尾部20 m，40 m，60 m，80 m，100 m处测量热尾流中心温度，如图3所示。可以看出，当网格数量大于516万时，尾流中心温度不再随网格数量增加而改变。因此，为节约计算时间，本文选择数量为516万的网格进行计算。

1.3 流动传热控制方程

质量方程

$\frac{{\partial \rho }}{{\partial t}} + \frac{\partial }{{\partial {x_i}}}\left( {\rho {u_i}} \right) = 0,$

(1)

动量方程

$\frac{\partial }{{\partial t}}\left( {\rho {u_i}} \right) + \frac{\partial }{{\partial {x_j}}}\left( {\rho {u_i}{u_j}} \right) = \frac{{\partial p}}{{\partial {x_i}}} + \frac{{\partial {\tau _{ij}}}}{{\partial {x_j}}} + \rho {g_i},$

(2)

能量方程

$\begin{split}&\frac{\partial }{{\partial t}}\left[ {\rho \left( {{c_p}T + \frac{{{u^2}}}{2}} \right)} \right] + \frac{\partial }{{\partial {x_j}}}\left[ {\rho {u_j}\left( {{c_p}T + \frac{{{u^2}}}{2}} \right)} \right] = \\ & - \frac{\partial }{{\partial {x_j}}}\left( {p{u_j}} \right) + \frac{\partial }{{\partial {x_j}}}\left( {k\frac{{\partial T}}{{\partial {x_j}}}} \right) + \frac{\partial }{{\partial {x_j}}}\left( {{u_j}{\tau _{ij}}} \right){\text{。}}\end{split}$

(3)

其中：

${\tau _{ij}} = \left[ {\mu \left( {\frac{{\partial {u_i}}}{{\partial {x_j}}} + \frac{{\partial {u_j}}}{{\partial {x_i}}}} \right)} \right] - \frac{2}{3}\mu \frac{{\partial {u_i}}}{{\partial {x_i}}}{\delta _{ij}}。$

(4)

第p相的体积分数满足方程：

$\frac{\partial }{{\partial t}}\left( {{\alpha _p}{\rho _p}} \right) + \frac{\partial }{{\partial {x_i}}}\left( {{\alpha _p}{\rho _p}{u_i}} \right) = 0.$

(5)

式中：ρ为密度，kg/m³；u_i（i=1，2，3）分别为x，y，z方向速度分量，m/s；p为流体静压，Pa；g_i为重力加速度分量，m/s²；c_p为定压比热容，J/（kg·K）；T为温度，K；k为导热系数，W/m·K；μ为粘性系数，Pa·s；α_p为第p相体积分数。

VOF法是目前求解分层流动、有自由表面流动等存在2种或多种互不相溶流体的复杂多相流动问题中最常用的方法。它在计算域内对互不相溶流体求解同一个动量方程，并追踪每种流体的体积分数来模拟多相流动。计算单元内各相体积分数之和为1，在单元内，若α_p=0，则单元中没有第p相流体；α_p=1，单元中充满第p相流体；0<α_p<1，单元中有第p相流体与其他相流体的界面。根据局部α_p值，计算域内每一控制容积被赋予适当的物性和变量值。湍流流动是非常复杂的流动，在计算湍流流动时需要附加湍流方程，本文研究均采用标准k-ε湍流模型。

1.4 物性参数及边界设置

参照海水物性参数随温度变化数据，拟合得计算所用相应物性参数关系式如下：

$\rho {{ = }}644.895 + 2.612T - 0.005{T^2},$

(6)

${c_{{p}}}{{ = }}6\ 545.302 - 15.377T + 0.025{T^2},$

(7)

$k{{ = }} - 0.666 + 0.007T - 8.699 \times {10^{ - 6}}{T^2},$

(8)

$\mu = 0.029\ 303{{ - }}0.000\ 132T + 4.109 \times {10^{ - 10}}{T^3}{\text{。}}$

(9)

潜艇排放冷却水与螺旋桨的搅动作用效果可视为“热射流”^[1]。针对目标潜艇，其热功率为220 MW，若热量全部通过冷却水散出，假定海水温度为284 K，冷却水温度为318 K，海水热容取4 000 J/（kg·K），则排水速度为0.52 m/s；假定潜艇航速为1 m/s，海水入口速度为0.1 m/s；为减弱壁面影响，将海水出口及两侧面定义为自由流；壁面法中将海面与海底均定义为壁面条件；VOF法中海面上空增加10 m空气层，海面为空气与海水交界面，空气层上表面给定4 m/s水平速度。

2 计算结果 2.1 壁面法热尾流浮升扩散规律

图4为纵截面上热尾流温度分布云图。以潜艇尾部排水口为原点，沿潜艇航行反方向，测量距潜艇尾部不同距离上热尾流中心区温度及其所在高度，可得尾流中心温度与浮升高度随其在水平方向扩散距离变化曲线如图5和图6所示。为方便分析，将图6中浮升高度曲线与水平方向（图中X轴方向）的夹角定义为热尾流的“浮升角”。显然，浮升角越大，热尾流浮升速率越高，当浮升角为90°时热尾流垂直向上浮升，此时浮升速率最高；当浮升角为0°时浮升速率最低为0。

观察图5和图6可得，在距离小于15 m区域内，由于热尾流刚从潜艇尾部排出，与周围环境存在较大的温差，会迅速与环境发生剧烈的能量、动量交换，尾流温度迅速降低；这阶段浮升角处于最大，浮升高度迅速增加；随着距离增加，尾流与海水间温差逐渐降低，热量交换逐渐减弱、浮升角逐渐减小，热尾流温度衰减、浮升速率逐渐降低；到达水面后，浮升角度减小为0°，热尾流不再浮升，在水面逐渐扩散、衰减，最终形成水滴形温度轮廓，如图7所示。据此，可以将热尾流水下发展过程大致分为：近尾流快速扩散、浮升区、过渡区、远尾流水面扩散区3部分。

图8分别为3~5 μm及8~12 μm波段下水面红外识别图像，可以看出，当温差大于0.06 K时，热尾流具有较为明显红外识别轮廓；在8~12 μm波段下，热尾流红外成像更加清晰，识别效果更佳。

2.2 VOF法热尾流浮升扩散规律

海面是大气与海洋相接形成的一个动态变化的气液交界面，它对大气与海洋之间热量、动量交换等具有重要的影响。为研究获得更接近于真实海况条件下的热尾流浮升规律及水面温度分布特性，本文采用VOF法对海面空气流动、气液交界面对热尾流的影响进行了研究。

图9为采用VOF方法获得的纵截面热尾流温度分布云图。图10为2种方法下尾流浮升高度随水平扩散距离变化曲线，图中给出了VOF法中热尾流浮升过程中3个区域划分。可以看出，在近尾流区2种方法热尾流浮升速率基本相等，但是由于VOF方法中不再将海面简化为壁面边界，而是随海水一起流动的气液交界面，受海面空气流动及交界面作用，海水流动加强，换热量增加，导致过渡区热尾流浮升角减小，浮升速率降低，从而使得浮升过渡段增长，浮升过程能量衰减增加；当浮升至水面后，受空气流动及交界面的影响，水面温度场出现明显抖动失稳现象，导致水面温度轮廓不再连续，如图11所示。

为方便分析，将图11中由热尾流引起的水面高温区域称为“温度异常区”。图12为温度异常区与环境海水温度之差大于某一值所占区域面积（A_T）的折算尺寸（将温度异常区转换为与其面积相等的正方形的边长长度（L_c），即L_c²=A_T）变化曲线。可以看出，由于受海面空气流动及海-气交界面的作用，热尾流浮升速率降低，浮升过程温度衰减增加，从而导致VOF法中浮升至水面的温差特征减弱，在壁面法中温差大于0.06 K的区域折算尺寸可达201 m，而在VOF法中减小了约55%，仅为91 m。可见，海面空气流动及海-气交界面对热尾流的浮升过程及表面温度分布特征具有重要的影响，VOF法对于有自由表面的多相流动具有较高的计算精度，仿真研究中可采用VOF法研究海面空气流动及交界面对热尾流的影响。

图13为VOF法在3~5 μm及8~12 μm波段下水面红外识别图，可以看出，热尾流依然具有较好的红外识别性。

2.3 不同排水方案热尾流浮升扩散规律

潜艇冷却水出口温度、出口速度对热尾流浮升、扩散过程及表面温度分布特征都具有重要影响。潜艇正常航行时，冷却水出口速度与温度应满足关系式：

${v_1} = \frac{Q}{{\rho {c_p}s({T_1} - {T_0})}},$

(10)

式中：T₁和T₀分别为冷却水出口温度与环境温度。

可以看出，相同热工况下，冷却水排出速度与温度成反比，增加冷却水出口温度，则出口速度会减小。当冷却水出口温度增加时，热尾流与环境温差增加，热尾流温度信号增强；当冷却水排出速度减小时，则单位时间内排出的热量减小，热尾流信号减弱。因而，理论上存在一组冷却水排出温度与速度的最佳值，可以使潜艇热尾流引起水面温度信号最弱、表面热红外特征最不明显。如表2所示，本文选取以下5组冷却水排放方案，研究冷却水排出速度与温度对热尾流浮升及水面温度分布特征的影响规律。

图14为不同排水方案下热尾流浮升高度随距离的变化曲线，可以看出，随着排水温度升高，热尾流浮升速率逐渐增加。这是由于，冷却水出口温度增加后，其与周围环境的温差增大，密度差增大，产生的浮升力增大，从而加快了热尾流浮升速率。

图15为不同排水方案下温度异常区的折算尺寸变化曲线。可以看出，随着排水温度升高，热尾流引起的水面温度异常区面积逐渐增加，排水温度由313 K增加至333 K，温差大于0.06 K区域折算尺寸由60.5 m增大至109.5，增加约81%。这是由于排水温度升高后，热尾流浮升速率增加，浮升过程能量衰减减小，从而导致水面的热尾流温度信号增强。

通过以上分析可得，排水温度变化对热尾流浮升过程及表面温度分布特性的影响要大于速度变化的影响，潜艇航行时可通过适当增加冷却水排出速度，降低排出温度，以减弱潜艇尾迹信号强度，增加自身隐身性能。

3 结　语

本文采用动网格技术对潜艇热尾流浮升过程及其水面温度分布与红外识别特征进行了仿真分析，研究表明：

1）潜艇热尾流浮升、扩散过程可分为：近尾流快速扩散、浮升区、过渡区、远尾流表面扩散区三部分；当温差大于0.06 K，热尾流具有较为明显红外识别轮廓。

2）海面空气流动及海-气交界面对热尾流浮升过程及其表面温度分布具有重要影响：受空气流动及交界面影响，热尾流浮升速率降低，水下漂浮距离增加约68%，表面温差大于0.06 K区域尺寸减小约55%。仿真研究中可采用VOF法研究海面空气流动及交界面的对热尾流的影响。

3）排水温度对热尾流浮升过程及其表面温度分布具有重要影响，排水温度增加20 K，表面温差大于0.06 K区域尺寸增大约81%。潜艇航行时可通过适当增加冷却水排出流量，降低排出温度，以减弱潜艇尾迹信号强度，增加自身隐身性能。