蔡君蕾, 姚天成, 刘 宏, 万立健 万 军, 樊 翔, 赵永生

(1. 上海交通大学 海洋工程国家重点实验室,上海 200240; 2. 自然资源部第二海洋研究所 极地深海技术研究院, 浙江 310012; 3. 上海交通大学 海洋装备研究院; 高技术船海装备创新研发中心,上海 200240; 4. 上海达华测绘科技有限公司,上海 200120;5. 上海船舶研究设计院,上海 201203)

海洋观测技术是海洋科技格局发展中最具活力的创新领域之一.中国作为海洋大国,对海洋环境的监测和水文资源要素的采集有迫切需求.传统的海洋环境监测通常以浮标、潜标和载人调查船为载体,存在难以适应复杂恶劣的海况、机动性能较差、综合成本较高等方面的天然局限,因此海上监测载体正在向无人化、智能化发展.无人艇[1]作为一种无人海洋智能运载平台,可在海洋中承担长时间、大范围、低成本的海洋监测任务,其中双体无人测量艇因其优良的横稳性、快速性及操纵性能和充足的甲板面积,正成为研究开发的热点方向.

在无人测量艇的设计过程中,航行阻力大小直接决定了测量艇的经济性能和工作性能,因此阻力性能的预报十分重要.目前,阻力预报方法主要有模型试验、理论计算和计算流体力学(CFD)方法.罗良[2]利用Holtrop经验公式法对KCS船模和实船分别进行分析研究,并开展物理试验进行对照,最终证明该方法能满足工程精度要求;孙源等[3]利用半理论半经验方法对某滑行艇进行阻力预报,并与物理试验结果进行比较,得出在高速滑行区间,阻力结果与试验值较接近,但纵倾角的计算值还存在一定误差;邓芳等[4]、杨显原等[5]运用细长体理论对双体无人船进行了阻力预报,与调研的参考文献试验结果进行验证对比,分析不同参数在减阻方面的影响;董文才等[6]比较分析了滑行艇阻力计算中常用的2类半经验半理论方法的适用范围,提出查洁法结合雷诺平均(RANS)方程的滑行艇阻力计算方法,与物理试验值进行比较,结果吻合较好.由此可知,物理试验是阻力预报最可靠的方法,但由于成本和试验设施因素,在船舶设计初期通常难以开展;而基于经验公式和势流理论的计算方法在阻力预报方面具有一定应用价值,但在精度方面仍有所欠缺.

CFD方法可以考虑流场的黏性影响,模拟更真实的实际流场情况,且相较于经验公式和势流理论方法,CFD方法能达到更高的精度,解决很多黏性流场带来的非线性问题,如对极端海浪、甲板上浪、波浪破碎等现象中船舶运动和载荷响应预报的模拟.不少学者也利用该方法进行阻力研究.倪崇本等[7]对阻力成因进行分析,提出用基于CFD方法的叠模法求解黏性阻力,用Euler方程求解兴波阻力来进行总体阻力预报;肖国权等[8]对小型单体无人船进行了CFD阻力分析,得到阻力关于不同参数的变化;邵峰等[9]利用FLUENT模拟了带有附体的小型无人船阻力预报,数据吻合良好;方静等[10]利用OpenFOAM对一艘超小型双体无人船进行CFD阻力分析,得到阻力性能较优的船体线型和非常规桨的布置方案;戴现令等[11]、梁严等[12]也分别利用FLUENT对三体无人船进行了阻力分析,来对片体型线、位置布置等参数进行阻力性能优化;赵核毓等[13]用CFD方法对深V滑行艇的多个船型进行阻力预报比较,为后续研究提供了数据参考;郭军等[14]提供了一种Savitsky方法结合CFD重叠网格技术的预报方法,提高了CFD方法对滑行艇阻力预报的计算效率和精度.综上所述,当前对小型多体无人船的阻力研究主要集中在中高航速区间的阻力预报,而对无人双体船从低航速到高航速区间和涉及较完整的航态变化的阻力性能研究还较为有限.

1 目标船总体方案

目标艇为上海交通大学船舶与海洋工程设计研究所自主研发设计的新型多功能无人双体测量艇.船体片体横剖面采用双折线设计,片体型线关于片体中心线不完全对称.该测量艇主要用于搭载多波束声学设备,也可以搭载其他探测设备进行海上作业.该船的外观效果图如图1所示,船体横剖面型线如图2所示,实船及船模的主要设计参数如表1所示.

表1 主要设计参数Tab.1 Main parameters

图1 目标船外观效果Fig.1 Appearance rendering of target ship

图2 目标船横剖面型线Fig.2 Cross-section of target ship

2 船模水池静水阻力试验

试验在上海交通大学多功能船模拖曳水池完成,水池全长 300 m,宽16 m,深7.5 m,并配备最大速度为10 m/s的高速拖车.

静水阻力试验在设计吃水状态下进行,试验过程中,船模通过连接在模型重心位置的前后共4根钢丝拖曳前行,钢丝与模型连接处的4个拉力传感器的合力即为模型总阻力.试验中开放纵倾和升沉两个自由度,通过倾角仪测量得到模型的纵倾角度,通过Qualisys光学测量系统记录升沉并换算至模型重心位置,同时测量记录水池温度与拖车速度.

图3 试验整体布置Fig.3 Overall layout of model test

3 目标船静水阻力性能预报

在数值计算方面,采用两种阻力预报方法进行比较.其中,基于统计学回归分析的Holtrop经验公式法依托Maxsurf平台[15]进行,主要进行阻力成分的分析比较;黏性计算流体力学方法则借助STAR-CCM+平台开展,用于全航速的阻力和航态模拟,最后将各个预报结果与水池物理试验结果进行对比和分析.

3.1 Holtrop经验公式法

Holtrop方法由荷兰 MARIN 水池在1978—1984 年发表,是根据船模试验和实船试航资料基于统计学方法分析而得到的实船阻力计算公式.其数据资料来自非系列排水型船舶,且涵盖弗劳德数Fr在0.55以上的高速船,可用于预测油轮、杂货船、渔船、拖船等不同船型的阻力,适应范围较广.阻力一般估算公式[16]为

Rt=Rf(1+k)+Rapp+Rw+Rb+Rtr+Ra

(1)

式中:Rf为摩擦阻力,根据ITTC1957公式[17]得到;1+k为船舶形状因子;Rapp为船舶附体阻力;Rw为兴波阻力;Rb为球鼻艏引起的附加阻力;Rtr为方艉引起的附加阻力;Ra为船模与实船相关的修正因子.Holtrop方法[18-19]给出了式(1)中每个参数的回归公式.

由于研究对象为双体船,故该方法中的形状因子可表示为1+βk,其中β为双体船片体间的黏性干扰因子.通过Molland公式[20]来定义形状因子,

(2)

3.2 黏性计算流体力学方法

黏流CFD方法选取的控制方程为RANS方程,数值离散方法为有限体积法,湍流模拟选取k-ε湍流.针对本次模拟,采用VOF方法对自由面进行捕捉,同时使用重叠网格技术和DFBI模型来模拟运动船体.

3.2.1控制方程 流场模拟基于RANS方程进行,在欧拉坐标系下,两相不可压缩流体运动的基本控制方程[21]为连续性方程和动量守恒方程,即N-S方程,分别为

(3)

(4)

为使上述方程封闭,需要引入湍流模拟(方程),根据文献[22]中的研究,选取了可实现的k-ε两层湍流模型,该模型在鲁棒性和准确性之间提供了良好的折衷,能更好地模拟船体流线弯曲程度较大的流动.

3.2.2基于重叠网格技术的DFBI模型 重叠网格是将复杂的流动区域分成多个几何边界比较简单的子区域,各子区域中的计算网格独立生成,彼此存在重叠、嵌套或覆盖的关系,流场信息通过插值在重叠交界面上进行匹配和耦合.在处理物体非定常运动数值仿真时,重叠网格可以解除物体与网格之间的拓扑结构约束,物体在计算域内自由运动,计算网格无需重新生成,将该方法应用于物体运动幅度较大的数值模拟中,能提高模拟结果的精度和稳定性.DFBI是一种通过计算刚体在流体中的力来推导刚体轨迹的模型,建模时可根据需要分别设置刚体6个自由度的约束,以进行针对性的自由度响应分析.结合重叠网格技术建立DFBI模型,可实现流体作用下刚体运动特性更真实的模拟.

3.2.3网格划分与边界条件 在STAR-CCM+平台中建立模型尺度的数值波浪水池模型,数值波浪水池长10 m,宽8 m,高度6 m.边界设置如图4所示,在距离压力出口2 m以内的区域设置消波区,即一个船长的阻尼消波区,保证达到消波效果又尽量不影响船舶运动时产生的兴波流场.

图4 数值波浪水池边界设置Fig.4 Boundary setting of numerical pool

网格划分采用切割体网格,主要在自由液面的兴波部分进行加密,计算域网格划分如图5所示.

图5 计算域网格划分Fig.5 Grid division of computational domain

表2 网格收敛性测试阻力结果(船模航速3 kn)Tab.2 Results of grid convergence study about resistance (model speed=3 kn)

由表2可知,粗糙网格得到的阻力结果为11.23%,误差较大,而中等网格和细化网格得到的结果与水池试验值的相对误差均在7%左右,在可接受范围.因此,综合考虑计算精度和计算效率,最终选择中等网格的划分方式.

4 阻力计算与试验对比分析

图6 静水阻力预报结果Fig.6 Results of hydrostatic resistance predictions

图7 阻力分量变化Fig.7 Resistance component of ship model

图8 CFD-物理试验纵倾值Fig.8 Trim of ship model (CFD-model tests)

图9 CFD兴波俯视图Fig.9 Vertical view of waveform around ship model by CFD

综上可知,Holtrop经验公式法主要适用于测量艇排水状态的阻力预报,对黏性流场的模拟较粗糙,与真实航行状态相比虽存在一些偏差,但计算效率高,在船舶设计初期可以快速预报测量艇在低速测量工况下的阻力性能;而黏流CFD方法能完整模拟预报测量艇的所有航行状态和阻力值,且能达到较高精度,但需花费较长时间,在32核数(CPU主频2.5 GHz)的计算机上进行单个工况的模拟所需时间约为8 h.

5 测量工况下的航态优化分析

测量艇航态是影响其阻力性能的关键因素,选取适当的初始航态能起到良好的减阻效果.一般在设计排水量下,通过适当调节压载分布来改变重心的纵向位置是改变航态的有效方法,以及通过对比目标艇不同重心纵向位置下的阻力性能进行航态优化分析.计算工况选取重点关注的实船工作航速 6 kn,设计了10个不同重心纵向位置的工况,如表3所示.

表3 航态优化分析工况Tab.3 Conditions for attitude optimization

工况5为初始设计重心纵向位置,在此基础上选取4个尾倾工况和5个首倾工况,得到的阻力和纵倾值关于重心纵向位置的变化如图10所示.

图10 航态优化阻力与纵倾值变化Fig.10 Resistance and trim of attitude optimization

由图10可知,通过改变重心的纵向位置,即主要改变船舶的纵倾航态,对静水阻力影响十分显著.在这10个工况中,重心纵向位置的移动范围为平衡位置处约±15%的船长,静水阻力最大值和最小值有35%的相对变化,其中重心的纵向位置在3.6 m和4.8 m处阻力达到了两次极小值,说明这是测量速度下两个比较经济的航态,这两个工况下船体周围波形如图11所示.

图11 两种航态船体周围波形Fig.11 Waveform around ship model of two navigation posture

由图11可知,在阻力值相似的条件下,重心纵向位置在3.6 m处时,船体产生的开尔文波波幅较小,船头波形也相对平稳,没有太多的波浪破碎现象;重心纵向位置在4.8 m时,船体产生的开尔文波幅较大,存在约6° 的首倾,即出现了埋首现象,船头伴有部分波浪破碎现象,整体稳定性不如前者.故经综合考虑,选定重心纵向位置在距离船尾3.6 m,即距离船尾45%的船长处,纵倾为1.4° 的尾倾为船模测量速度下的经济航态,此时船体所受阻力较小,开尔文波形和纵向航态均趋于平稳,适合需要稳定环境下的探测设备进行工作.

6 测量传感器布置位置优化

目标艇用于搭载传感器进行海洋观测作业,传感器安装在测量艇的两个片体中间开槽区域,具体安装位置根据不同工况决定.以多波束声学设备为例,该探测设备的安装位置需保证传感器所处介质的统一,避免传感器出入水,故需在高液位区域安装传感器,保证其工作中始终处于自由液面以下.为此,在选取经济航态基础上,用CFD方法监测实船工作航速下的静水阻力片体间波形,5和6 kn航速下片体间波形的左视图和俯视图分别如图12和图13所示.

图12 实船5 kn航速下的片体间波形Fig.12 Waveform between hulls at 5 kn

图13 实船6 kn航速下的片体间波形Fig.13 Waveform between hulls at 6 kn

由图可知,航速5 kn时,片体间高液位峰值在距船尾0.47个船长处;航速6 kn时,片体间高液位峰值在距船尾0.39个船长处.综合考虑得出,在工作航速时,距船尾约0.3～0.5个船长的位置存在一个相对稳定的高液面区域,如果把多波束声学探测设备放在此处,可有效减少设备出水风险;同时,高液位区可远离自由面,以避免自由面上的气泡和乱流影响.综上所述,可以确定该工作航速下传感器的推荐安装位置,根据不同的工况和设备需求也可通过同样的方法调整设备的安装位置.

7 结论

综合利用经验公式、CFD方法和水池物理模型试验等方法对自主研发的双体无人测量艇进行静水阻力性能研究,同时分析不同航态下的阻力变化和工作航速下片体间流场形态,得到主要结论如下:

(2) 本双体测量船在重心纵向位置从平衡位置移动±15%船长的范围内,静水阻力有35%的相对变化,说明纵倾值对阻力性能的影响十分显著.在测量航速下得出本测量艇较为经济的航态为1.4° 的尾倾,此时重心纵向位置在距离船尾45%船长处,由此可知与本船型相近的双体船的推荐初始航态为1.5° 左右的尾倾,重心的推荐纵向位置布置点位于船体中部偏后处.

(3) 本船工作航速下,片体间距离船尾约0.3～0.5个船长的位置存在一个相对稳定的高液位区域,故多波束声学探测设备的推荐布置位置为片体间船体中后部处,以保证多波束设备的探测精度和稳定性.