螺旋桨诱导的船体表面脉动压力数值预报*

阳 涛 林永水 吴卫国 秦江涛 贺 伟 秦升杰

(武汉理工大学交通学院1) 武汉 430063) (武汉理工大学理学院2) 武汉 430070) (国家深海基地管理中心3) 青岛 266237)

螺旋桨诱导的船体表面脉动压力数值预报*

阳 涛1)林永水1,2)吴卫国1)秦江涛1)贺 伟1)秦升杰3)

(武汉理工大学交通学院1)武汉 430063) (武汉理工大学理学院2)武汉 430070) (国家深海基地管理中心3)青岛 266237)

基于RANS方法，采用滑移网格技术建立“船+桨”计算模型对双尾客船的螺旋桨诱导的脉动压力进行了数值预报.对比不同航速下船模阻力试验值与计算值，发现二者吻合良好，计算精度较高，说明了计算模型的可靠性；重点分析了船底板脉动压力沿船长方向及船宽方向的分布规律，发现最大脉动压力分别发生在桨盘面前方0.1D处以及桨轴心线右侧0.05D处，且脉动压力沿船长方向衰减速度较船宽方向更快，而在船宽方向的作用范围较船长方向更广；分别探讨了不同梢隙比、不同转速、不同航速下脉动压力的变化规律.研究表明，脉动压力一阶幅值是该船后螺旋桨诱导脉动压力的主要成分，增加梢隙比、减小螺旋桨转速，以及增加航速均会降低螺旋桨诱导的脉动压力.

螺旋桨；脉动压力；数值预报；船体阻力；网格划分

0 引 言

自20世纪70年代以来，随着船舶不断向大型化和高速化发展，船舶主尺度、载重量，以及主机功率均在不断加大，然而受到尾型、吃水，以及水动力性能等条件的限制，螺旋桨直径及个数增加空间很小，这使得桨叶负载与转速不断增加，伴随着船尾不均匀流场的存在，由螺旋桨所诱导的激振力持续增大，最终导致船尾螺旋桨上方振动骤增，船尾有害振动事件频发，据统计，船舶上出现的有害振动，80%的案例是由螺旋桨激振力所应起的[1].螺旋桨在船后不均匀的伴流场中运转，随着桨叶负荷周期性改变，必将对船体产生激振力.根据作用力的传递方式不同，螺旋桨激振力可分为轴承力与表面力：轴承力是指以叶频为基频的脉动力和力矩通过轴系传递到船体的作用力，表面力即由螺旋桨诱导的船尾压力场经流体传递至船体表面的激振力.文献[2]中的实测结果均表明螺旋桨诱导的脉动压力诱发船尾有害振动的最主要因素，这也是引起国内外学者研究的重要原因所在.

目前，国内外螺旋桨脉动压力的计算方法主要有理论计算、经验公式、试验方法和CFD数值计算.理论计算主要基于流体力学理论，通常需要采取一些假设该条件，主要方法包括势流理论、升力面法、面元法，以及计算流体力学方法等.经验公式主要包括Holden法、Johnsson法、谷口中-高桥肇方法等[3]，各公式通常只适用于特定尾型的船舶，其精度较低，有一定的局限性.试验方法需要在空泡水筒或者减压拖曳水池中进行，对试验硬件设施要求较高，国内主要以熊鹰、伍锐、季斌等学者开展了较多研究.CFD数值计算具有低成本、高效快速等优势，随着近年来计算机软硬件技术的不断发展，CFD数值方法在计算精度与稳定性上有了较大的提升，第26届国际拖曳水池会议中CFD委员会明确指出CFD计算精度可以满足工程应用要求.殷玉梅[4]基于双向流固耦合，选用RNGk-ε模型对船桨流体模型进行脉动压力预报，并在德国汉堡船模试验水池(HSVA)开展模型试验，对比仿真、实验结果及经验公式结果，发现CFD方法与实验结果较为接近，其精度要远高于经验公式.丁科等[5]以E779A桨为研究对象，采用RANS方法结合VOF方法考虑自由液面，通过进行“船+桨”模型的数值自航试验，通过瞬态计算得到螺旋桨激振力的时域值，得到了船体表面脉动压力的大小.Sato等[6]应用CFD数值仿真对10个不同的螺旋桨进行了脉动压力的预报，并将预报结果与模型试验值相比较，结果显示，CFD方法能够较为准确的应用于螺旋桨片空泡脉动压力预报，同时对于脉动压力一阶幅值的预报已有较高的精度.

目前对于脉动压力的研究往往只针对较少测点位置进行分析，鲜有对其在船体表面的分布规律开展深入研究，较少的测点数据不足以为船舶振动计算及舱室噪声预报提供激励数据.基于此，文中通过在船底板布置丰富测点，对螺旋桨诱导的脉动压力进行了预报，对其在船底表面的分布规律进行了重点分析，同时探讨了脉动压力在不同梢距、转速以及航速下的变化规律.

1 计算模型与阻力验证

1.1 计算对象

本研究对象为一艘498客位现代客船，其为双机双桨单体船，主要参数见表1～2.CFD计算模型的缩尺比为1∶10.

表2 螺旋桨主要参数

1.2 计算域的建立

计算域的建立既要保证计算可靠性，也要考虑硬件的承载能力.文中研究对象为双机双桨的左右舷对称船舶，左右双桨之间的距离为3.6 m，经计算发现双桨均对彼桨上方D×D区域的流场干扰较小，而螺旋桨上方D×D区域为本文的主要研究区域，因此不考虑双桨干扰因素，只针对左半舷船舶与螺旋桨进行分析，以提高计算效率.

计算域整体为一个长5Lwl，半径为1Lwl的1/4圆柱体流域，其相较于长方体流域可以减少网格数量.其中，为保证来流的均匀性，设置流域入口至船首距离1Lwl；为尽量减小壁面效应的影响，设置1/4圆柱计算域半径1Lwl；为保证桨后尾流的充分发展，设置流域出口至船尾距离3Lwl.计算域整体见图1.

图1 计算域

为了后续网格划分时进行网格质量与数量控制的方便，整个计算域划分为流体大域、船体流域以及螺旋桨流域，见图2.其中，船体区域其宽与高均为2倍船宽，前后边界分别距离船体首尾距离均为1倍船宽；螺旋桨区域为一个与螺旋桨同轴的圆柱流域，将螺旋桨包裹其中，其直径为1.2D.

图2 船体流域与螺旋桨流域

1.3 网格划分

为更好地拟合船体与螺旋桨曲面、减小网格数量，船体及螺旋桨区域采用非结构化网格，流体大域采用结构化网格.对于船体流域，主船体表面网格尺寸设置为0.001 5Lwl，因艏部球首以及艉部表面曲率变化较大，网格划分应进行适当的加密；对于螺旋桨区域，桨叶表面网格尺寸设置为0.01D，而导边、随边、叶梢及叶根处是网格划分时应关注的重点，应进行适当加密.由于壁面网格尺寸的选取对于湍流模型的计算结果有着较大的影响，因此需要对船体表面近壁网格进行棱柱体网格划分，经多次试算后的Y+值对比，最终确定设置4层棱柱网格，其中第一层网格厚度设置为0.000 5Lwl，而后以1∶1的比例增长，以确保湍流模型最终的壁面Y+值在30～300之间.网格质量(quality)、纵横比(aspect ratio)以及最小角度(min angle)是进行网格划分时的主要监控对象，通过对网格进行多次光顺优化，最后使得以上3个指标分别控制在0.35°,0.25°,18°以上，总网格数目控制在293万个.

图3 船尾与螺旋桨网格

1.4 边界条件与计算设置

入口设置为速度入口(velocity-inlet)，V=2.44 m/s(对应实船航速15 kn)，出口设置为压力出口(pressure-outlet)，自由液面与对称面设置均为对称表面(symmetry)，流体壁面设置为移动壁面(moving-wall)，船体表面、桨轴及支架、螺旋桨表面均设置为静止壁面(stationary-wall).计算模型采用瞬态计算，设置螺旋桨流域网格为滑移网格，绕桨轴额定转速为1 600 r/min.螺旋桨流域与船体流域之间，以及船体流域与流体大域之间均采用Interface连接，各流域之间通过Interface插值进行信息传递[7-11].湍流模型采用RNGk-ε模型/标准壁面函数，压力速度耦合方式设置为PISO，计算时间步长为0.000 469 s(每步旋转4.5°).

1.5 船体阻力验证

在船桨模型进行脉动压力预报的过程中，可对船体表面阻力进行监测(总阻力的一半).因模型采用瞬态计算，因此船体阻力不可能像稳态计算达到绝对收敛，只能达到相对收敛.以实船航速Vs=15 kn(Vm=2.44 m/s)为例，见图4a)，瞬态计算的模型阻力(已由半船换算至整船)在75.5～75.9 N之间变化，其平均阻力为75.785 N，将计算值与试验值进行对比，发现二者误差为2.499%.与此同时，后续开展对不同航速下的螺旋桨脉动压力预报，相应得到各航速下船体受到的阻力.对比不同航速下船体阻力计算值与实验值，见图4b)，在各航速下(为方便对比，以实船航速为横坐标)，计算模型对于阻力的预报最大误差不超过4.755%，证明计算模型有较高的可靠性.

图4 船体阻力性能

2 螺旋桨脉动压力计算与结果分析

2.1 螺旋桨诱导的船体表面脉动压力分析

图5为船体表面脉动压力监测点设置.为分析螺旋桨在尾部船底表面所诱导的脉动压力分布规律，在螺旋桨正上方船底表面选取D×D区域进行脉动压力的分析.以螺旋桨正上方船底板位置为原点建立坐标系，设置船长、船宽、型深方向分别为X,Y,Z轴，并在D×D的船底板上布置有纵横8×7共56个的脉动压力监测点.其中Px4处Y向7个监测点均为探讨最大脉动压力而设置，其X向位置为X=0.1D，其他各监测点均以0.166 7D等距排布.

图5 船体表面脉动压力监测点设置

结合Matlab对Journal文件对监测点的坐标信息以及监测信息进行编写，并通过Fluent读取Journal文件进行监测点坐标的批量输入，以及批量监测设置.经过近72 h的计算，观察到各监测点压力变化出现明显的周期性，计算收敛，船体表面Y+值控制在187，螺旋桨叶背与叶面Y+值分别控制在171和154，近壁网格精度较高，见图6.

图6 船体与螺旋桨表面Y+值分布

选取计算收敛后螺旋桨旋转1周各监测点的时域压力数据进行分析，见图7.图7a)为桨轴心线上方船长方向8个监测点(P0，Px1～Px7)压力时域变化曲线，压力变化幅度最大的监测点出现在桨盘面前方0.1D处的Px4，由该点向船长方向前后延伸，各点的压力变化幅度逐渐减小；同时从桨后Px3处至桨前Px6，其平均压力逐渐降低，Px7处开始提升，这是因为螺旋桨的旋转在桨上方形成了一个低压区.图7b)为桨盘面正上方船宽方向7个监测点(P0，Py1～Py6)的压力时域变化曲线，其中P0点位于螺旋桨正上方，也是离叶梢最近的监测点，很明显能看出从螺旋桨正上方P0点向船宽方向左右延伸，各点的压力变化幅度逐渐减小，同时各监测点的静压变化，出现了由左至右的相位延时，这是因为螺旋桨旋转方向为右旋，叶梢靠近各监测点的时间是从左至右依次的.

图7 船体表面监测点压力时域曲线

选取在收敛后螺旋桨旋转12周的时域数据进行FFT处理，得到各监测点脉动压力频域变化曲线，见图8.由图8可知，各监测点脉动压力均在叶频(133 Hz)与倍叶频(266 Hz)处出现峰值，而在其他频率其幅值均接近于0 Pa(为显示方便，400 Hz以后曲线未予显示)，同时叶频处的峰值远大于倍叶频处的峰值，是螺旋桨诱导脉动压力的主要成分.

图8 船体表面监测点脉动压力曲线

图9为船长方向(Y=0处)与船宽方向(X=0处)系列监测点的脉动压力一阶幅值分布曲线.在船宽方向，为更准确获取最大脉动压力幅值及其发生的位置，在P0与Py1之间增设2个监测点.由图8～9可知，在船宽方向，对于沿轴心线左右对称的位置，右侧位置的脉动压力均要大于左侧，而最大脉动压力发生在桨轴心线右侧约0.05D(Y=-0.05D)的位置；在船长方向，最大脉动压力的位置并不出现在桨盘面正上方，而是在桨盘面前方约0.1D处，这与文献[1]中多条船舶的实测规律是一致的,这也间接验证了文中分析方法的可靠性.

图9 脉动压力分布曲线

为更全面的获得船底板上的脉动压力分布规律，本文对所设56个监测点经FFT处理所得的脉动压力一阶幅值进行插值拟合，得到螺旋桨上方底板D×D区域脉动压力分布云图以及三维视图，见图10.由图10可知，螺旋桨在船底板上诱导的最大脉动压力约分布于以Y≈-0.05D，X≈0.1D为中心的椭圆区域内；在船长方向，以X=0.1D为中心向两端延伸，脉动压力快速下降，至X=-0.2D以及X=0.5D处，脉动压力已由最大377 Pa下降至50 Pa以下；在船宽方向，脉动压力以Y=-0.05D为中心向两侧下降，其下降速度较船长方向慢，而其作用范围较船长方向更广.

图10 螺旋桨上方D×D船底板脉动压力分布

2.2 桨盘面不同梢距脉动压力计算分析

充分了解桨盘面不同梢距处脉动压力的分布规律，对于船尾的尾型设计有着重要的指导作用.由图11可知，在桨盘面正上方Zt/D=0.04，0.08，0.12，0.145(船底表面)处设置4×7个监测点(Zt为桨叶正上方监测点与叶梢之间距离).

图11 桨盘面不同梢距监测点布置

图12为螺旋桨正上方4个监测点的压力时域变化曲线，经FFT处理可得相应监测点脉动压力频域曲线见图13，图14为不同梢距处脉动压力一阶幅值沿船宽的分布曲线.结合图13～14可知，随着梢距的增大，脉动压力一阶幅值逐渐减小，其下降的趋势逐渐减缓；同时，梢距越小，脉动压力的二、三阶幅值体现的越明显；梢距变化对螺旋桨正上方位置的脉动压力幅值影响较大，越靠近左右两侧，其影响逐渐减小，这与文献[8]中对于某桨的试验结果分析的规律是一致的.

图12 螺旋桨正上方监测点压力时域曲线

图13 螺旋桨正上方监测点脉动压力曲线

图14 船宽方向不同梢距脉动压力分布曲线

2.3 不同转速螺旋桨脉动压力计算分析

为探讨螺旋桨转速对其诱导的脉动压力的影响规律，在保持航速不变的情况下，本文选择了转速n=400，1 450，1 500，1 550，1 600 r/min 5个工况进行脉动压力分析.选择桨轴对称线上X=0.1D处的监测点Px4作为分析对象，各转速下Px4的压力时域变化曲线见图15，经FFT处理可得到其脉动压力一阶幅值，其随转速变化曲线见图16.保持航速不变，随着转速的增加，螺旋桨诱导的脉动压力一阶幅值大致呈线性增长，同时其二阶及更高阶次增幅可忽略不计.

图15 不同转速Px4压力时域曲线

图16 不同转速Px4脉动压力曲线

2.4 不同航速螺旋桨脉动压力计算分析

图17 不同航速Px4压力时域曲线

图18 不同航速Px4脉动压力曲线

3 结 论

1) 螺旋桨在其上方船底板D×D区域诱导的脉动压力：其大小主要由其一阶幅值决定，二阶及更高阶次幅值已远小于一阶幅值；在船宽方向沿轴心线左右对称的位置，右侧位置的脉动压力均要大于左侧，最大脉动压力出现在桨轴心线右侧约0.05D的位置；船长方向的最大脉动压力出现在桨盘面前方约0.1D处；脉动压力在船长方向的衰减速度较船宽方向更快，而在船宽方向的作用范围较船长方向更广.

2) 螺旋桨桨盘面内，随着梢距的增加，其诱导的脉动压力一阶幅值逐渐下降，下降的趋势逐渐减缓；同时梢距变化对螺旋桨正上方位置的脉动压力影响较大，越靠近左右两侧，其影响逐渐减小.

3) 保持航速不变时，随着螺旋桨转速的增加，脉动压力一阶幅值大致呈线性增加；保持转速不变时，随着航速增加，螺旋桨诱导的脉动压力逐渐减小，其减小的趋势逐渐放缓.

4) 研究结果可以为船尾结构振动控制提供指导，同时也可作为螺旋桨脉动压力引起的船体振动以及舱室噪声预报时的输入激励.

[1]何友声,王国强.螺旋桨激振力[M].上海:上海交通大学出版社,1987．

[2]林永水,吴卫国,翁长俭.20 000 t近海散货船有害振动诊断及治理[J].武汉理工大学学报(交通科学与工程版),2015(2):292-296.

[3]王蒙蒙.螺旋桨诱导的船体表面脉动压力预报方法研究[D].大连：大连理工大学,2006.

[4]殷玉梅.船体振动固有特性及螺旋桨脉动压力预报方法[D].大连：大连理工大学,2015.

[5]丁科,王永生,魏应三.船尾螺旋桨激振力的数值预报[C].第十一届全国水动力学学术会议暨第二十四届全国水动力学研讨会并周培源诞辰110周年纪念大会文集:下册，北京：海洋出版社，2012.

[6]SATO K, OHSHIMA A, EGASHIRA H, et al. Numerical prediction of cavitation and pressure fluctuation around marine propeller[J]. Report of Ship Research Institute,2009(2):167-172.

[7]李亮,王超,孙盛夏,等.船桨舵一体的螺旋桨激振力数值预报分析[J].武汉理工大学学报(交通科学与工程版),2015(4):768-772.

[8]KUROBE Y, UEDA T. On pressure fluctuations induced by cavitating propellers: basic experiments and examples of reduction method of pressure fluctuations[J]. Report of Ship Research Institute,1982,19:49-72.

[9]傅慧萍.船桨整体及螺旋桨诱导的船体表面脉动压力计算[J].哈尔滨工程大学学报,2009,30(7):728-734.

[10]PAIK K J, PARK H G, SEO J. RANS simulation of cavitation and hull pressure fluctuation for marine propeller operating behind-hull condition[J]. International Journal of Naval Architecture & Ocean Engineering,2013,5(4):502-512.

[11]JI B, LUO X, PENG X, et al. Numerical analysis of cavitation evolution and excited pressure fluctuation around a propeller in non-uniform wake[J]. International Journal of Multiphase Flow,2012,43(43):13-21.

Numerical Prediction of Pressure Fluctuation Included by Propeller on Hull Surface

YANG Tao1)LIN Yongshui1,2)WU Weiguo1)QIN Jiangtao1)HE Wei1)QIN Shengjie3)

(SchoolofTransportation,WuhanUniversityofTechnology,Wuhan430063,China)1)(SchoolofScience,WuhanUniversityofTechnology,Wuhan430070,China)2)(NationalDeepSeaCenter,Qingdao266237,China)3)

By using Reynolds-Averaged Navier-Stokes (RANS) equations, a ship-propeller model based on sliding mesh method is performed to numerically predict the propeller-induced pressure fluctuation for a double-stern passenger ship. It is found that hull resistances are predicted fairly well compared to the experimental results under different speeds, which validates the calculation models. The distribution of pressure fluctuation on hull bottom is analyzed in detail. Results show that the maximum of pressure fluctuation occurs at 0.1Din front of propeller disk in the direction of length, and 0.05Dto the right of propeller shaft axis in the direction of width. The pressure fluctuation decreases faster in the direction of length, while it presents a wider influencing scope in the direction of width. The variation laws of pressure fluctuation are discussed under series of tip clearance ratios, rotational velocities and ship speeds, respectively. The results show that the first-order amplitude of pressure fluctuation is dominant. In addition, the increase of tip clearance ratio and speed of the ship and the reduction of rotational velocity of propeller can lower the induced pressure fluctuations.

propeller; pressure fluctuation; numerical prediction; hull resistance; grid generation

2016-11-18

*中央高校基本科研业务费专项资金资助(173114003)

U661.1

10.3963/j.issn.2095-3844.2017.01.019

阳涛(1991—)：男，硕士生，主要研究领域为船舶结构振动与噪声预报及控制