非线性水动力导数的数值计算与研究*

赵小仨 徐海祥

(高性能船舶技术教育部重点实验室1) 武汉 430063) (武汉理工大学交通学院2) 武汉 430063)

非线性水动力导数的数值计算与研究*

赵小仨1,2)徐海祥1,2)

(高性能船舶技术教育部重点实验室1)武汉 430063) (武汉理工大学交通学院2)武汉 430063)

针对船舶的非线性运动难以界定和非线性运动难以预报的问题，以供应船为研究对象，采用CFD商用软件FLUENT，结合动网格技术对大振幅平面运动机构试验进行数值模拟，通过对比不同工况的流场压力云图，分析得出供应船水动力达到非线性的振幅范围.设计供应船非线性运动的试验方案，分别模拟不同频率时的大漂角斜航运动及大振幅的纯纵荡、纯横荡、纯首摇、组合运动，拟合得到接近零频率的非线性水动力导数.

非线性水动力；大振幅PMM试验；数值计算；供应船

0 引 言

船舶操纵性与船舶航行安全紧密相关，是重要的水动力性能之一.近些年，国际海事组织(international maritime organization, IMO)前后颁布了A.751(18)和MSC.137(76)号决议，针对船舶操纵性的问题提出了明确的要求，并建议各国政府机构按要求执行.SIMMAN 2008和SIMMAN 2014的研讨会，评估了CFD(computational fluid dynamics, CFD)方法预报船舶操纵性的能力.第25届ITTC(international towing tank conference, ITTC)操纵会议对现有的船舶操纵性预报方法做了总结.总之，船舶操纵性能越来越受到造船界的重视[1-3].

水动力导数对船舶操纵性的预报至关重要.目前，通过平面运动机构试验(planar motion mechanism test, PMM)确定船舶水动力导数是最可靠的方法之一.从SIMMAN 2008发布了针对三个标准船模进行的PMM试验的实验数据以来，国内外学者开始对CFD模拟PMM试验进行探究.Turnock等[4-12]用CFD软件模拟小振幅PMM试验，求取线性水动力导数；Toxopeus等[13-16]建立非线性水动力模型，模拟大振幅PMM试验，求取非线性水动力导数.

虽然许多学者对数值模拟PMM试验做了大量研究工作，但是迄今没有一个定量的标准来判断船模的运动是否达到非线性范畴，相关文章也较少.评判船舶的运动是否达到非线性，不仅取决于运动参数，还与船型等因素有关.文中将以供应船为研究对象[17]，通过数值模拟船模不同运动幅值的PMM试验，分析出供应船水动力达到非线性的运动幅值范围.在此基础上，设计试验工况，计算零频率附近的非线性水动力导数.

1 数学模型

研究船舶在大振幅下的操纵运动，用线性水动力模型很难准确的表达船舶所受到的水动力，为了更准确的描述船舶的运动，须考虑运动状态的非线性项[18].根据经验，在粘性类流体动力和力矩的泰勒级数展开式中保留至三阶项，对描述船舶在常速域中的运动已足够精确.

1) 流体惯性力(矩)XI，YI，MI

(1)

式中：u，v，r分别为船舶的纵向、横向速度与转首角速度.

2) 粘性类流体动力(矩)Xh，Yh，Mh

(2)

3) 水动力导数的无因次化 采用MMG模型系统建议的以Ld为参考面积进行无因次化.

2 平面运动机构试验

PMM试验是约束模试验的一种.是通过测量船模所受到的水动力和力矩，求得计算船舶操纵运动所需的各种加速度导数、速度导数和耦合导数.PMM试验有小振幅和大振幅之分，前者只能测定线性水动力导数，后者可以测定非线性水动力导数.本文主要模拟PMM试验的以下几种运动形式：定漂角斜航；纯横荡运动；纯纵荡运动；纯首摇运动；组合运动.

3 数值计算方法

3.1 计算模型

计算模型为一艘75 m供应船，缩尺比为1∶20.船模几何参数见表1.三维模型见图1.

表1 供应船模参数

图1 供应船模型

3.2 坐标系

坐标系见图2，O点位于船舯；X轴指向船首为正；Y轴指向左舷为正；Z轴正向依据右手定则确定.

图2 随船坐标系

3.3 划分计算域

计算域尺寸船首上游取1.5倍船长，船尾下游取3倍船长，船两侧取2倍船长，水深方向取8.3倍吃水，见图3.

图3 供应船的计算域

3.4 划分网格以及验证收敛性

为了保证网格质量，采用分块全结构化网格，并在首尾部以及呆木处进行适当加密.

以表1中的船模为研究对象，采用3种不同数量的网格模拟相同工况下的纯首摇运动，进行网格收敛性的验证.网格数分别为100万，200万和300万.图4给出了不同网格数计算得到的纵向力X、横向力Y和转首力矩M在一个周期内的曲线.由图4可知，网格数量从100万增至200万，计算结果变化明显；当网格数量从200万增至300万时，计算结果几乎不变.因此选择250万左右的网格，既能保证收敛性，又能节省计算资源和计算时间.

图4 不同网格数的计算结果对比

3.5 数值计算方法

1) 边界条件 ①入流与出流边界条件，入流面设为速度入口；出流面设为自由流出口，权重为1；②船体表面，在船体表面施加无滑移壁面条件；③自由面，考虑到供应船舶在定位过程中时低速航行，忽略自由面兴波的影响，将自由面设为对称面.

2) 定义动网格(dynamic mesh) 编写纯横荡与纯首摇运动的用户自定义函数(UDF)，并使之与FLUENT相关联；网格更新方法(mesh methods)选择网格光顺方法(smoothing)和动态层方法(layering).

3) 离散格式和求解算法 非稳态流动；压强插值格式选用标准格式；空间离散采用二阶迎风差分格式；时间积分方案采用一阶隐式；求解算法采用基于速度-压力耦合的SIMPLEC算法.

4 试验设计及计算结果分析

4.1 供应船非线性模型的振幅范围

计算工况见表2～3.其中，v′=va/Uc和r′=ra·L/Uc分别为纯横荡和纯首摇无因次化的速度幅值.

表2 纯横荡运动计算工况

表3 纯首摇运动计算工况

图5～6给出了不同运动幅值时的纯横荡与纯首摇运动所受横向力Y在一个周期的变化曲线.图7～8给出了纯横荡与纯首摇运动横向力Y高阶量的实际值与三阶拟合值的对比.从图中可以得出以下结论.

1) 随着横荡幅值和首摇角幅值的增大，横向力幅值明显增大，但横向力的相位几乎没有发生变化.

2) 当振荡幅度超过某一幅值，横向力已经不再是一阶正(余)弦的形式，即运动幅值越大，非线性表现的越强.

3) 由图6可知，当初始首向角为8°时，纯首摇运动表现为一阶正余弦的形式；当初始首向角为18°时，纯首摇运动已经不再是一阶正(余)弦的形式.

4) 由图7～8可知，当横荡运动的横向幅值y超过0.9 m(运动过程中的最大漂角约为31°)，当首摇运动的初始首向角超过30°时，横向水动力高阶量的实际值与三阶拟合值间的误差越来越大.

图5 75 m供应船纯横荡运动横向力Y随横荡幅值的变化

图6 75 m供应船纯首摇运动横向力Y随初始首向角的变化

图7 纯横荡运动Y方向高阶量的拟合值与实际值的对比

图8 纯首摇运动Y方向高阶量的拟合值与实际值的对比

因此，对于75 m供应船模，其非线性水动力模型的振幅范围约为漂角β大于10°小于30°.

4.2 流场分析

图9～10分别给出了75 m供应船首摇角幅值为8°和25°的纯首摇运动一个周期内的船体表面压力云图.压力变化的总体趋势是：纵向来看，从船首向船中压力逐渐递减，从船中向船尾压力逐渐递增；横向来看，随着首向角转动位置与方向的不同，最大压力的分布左右交替出现.通过对比图9～10可知，首摇幅值为25°时的压差比首摇幅值为8°的压差要大得多.

图9 75 m供应船首摇幅值为8°的船体表面压力分布图

图10 75 m供应船首摇幅值为25°的船体表面压力分布图

由图9可知，前半周期船体表面的压力分布情况与后半周期是关于船体中纵剖面对称的.图9b)表示首摇运动3/8周期时的压力分布图，它的左舷压力分布恰好与图9c)中右舷压力分布相同.出现该种现象的主要原因是：纯首摇运动的前半周期与后半周期的转首角速度大小相同方向相反，而且当首摇幅值为8°时，属于小振幅振动，符合线性假定，因此由运动输出的力在前半周期与后半周期是左右舷对称的.

由图10可知，前半周期船体表面的压力分布情况与后半周期已不再是关于船体中纵剖面对称.其主要原因是：当首摇幅值为25°时，属于大振幅振动，运动的输入与输出已不是线性关系.

4.3 非线性水动力导数的计算结果

以75米供应船为研究对象，根据4.1节所得结论，确定非线性水动力的计算工况.进而通过模拟直航、大漂角斜航,以及频率为0.2～0.4 rad/s的大振幅纯纵荡、纯横荡、纯首摇和组合运动，计算各运动模态下的水动力和水动力矩，然后求取不同频率时非线性模型的水动力导数，最终得到75 m供应船在零频率附近的各水动力导数[19].

1) 数据处理方案 以纯横荡运动为例，介绍本文的数据处理方案.其他运动形式的数据处理方法类似.根据纯横荡的运动规律，非线性水动力数学模型可变形为：

(3)

分别计算不同振动频率时船体所受的水动力X，Y，N，再用最小二乘法把水动力按上式拟合，分离出系数Ys1，Yc1，Yc3，Ms1，Mc1，Mc3，进而得到相应的水动力导数.

2) 直航运动 通过模拟0.4m/s～0.8m/s的直航运动，可求取纵向速度导数Xu，Xuu，Xuuu.图11给出了不同航速时纵向力的变化曲线.用最小二乘法对图11中的曲线进行拟合，再根据纵向水动力的三阶非线性表达式，可以直接得出3个纵向水动力导数，结果见表4.

图11 直航运动纵向力的拟合曲线

无因次量X'uX'uuX'uuu数值×10-3-0.42-2.370.21

3) 斜航运动 通过模拟10°～30°不同漂角下的斜航运动，可求取速度导数Xvv，Yv，Yvvv，Nv，Nvvv.图12给出了不同漂角时纵向力、横向力和力矩的变化曲线.与直航运动求取水动力导数的方法相同，得出斜航运动对应的五个水动力导数，结果见表5.

图12 斜航运动中纵向力、横向力及首摇力矩的拟合曲线

无因次量X'vvY'vM'vY'vvvM'vvv数值×10-20.50-2.95-1.63-7.820.27

表6 纯纵荡计算结果

表7 纯横荡计算结果

7) 组合运动 通过模拟不同频率下首摇幅值为25°漂角为12°的组合运动，可求取水动力导数Xvr，Xrr，Yvvr，Yvrr，Mvvr，Mvrr.计算结果见表10～12.

表8 纯首摇横向力计算结果

表9 纯首摇转首力矩计算结果

表10 组合运动纵向力计算结果

表11 组合运动横向力计算结果

表12 组合运动转首力矩计算结果

5 结 束 语

文中以供应船为研究对象，在FLUENT中数值模拟供应船模在不同运动幅值的PMM试验，进而对比各流场压力云图，分析得出供应船水动力达到非线性的振幅范围约为漂角β大于10°小于30°.最后，设计供应船非线性运动的试验方案，分别模拟不同频率时的大漂角斜航运动以及大振幅的纯纵荡、纯横荡、纯首摇、组合运动，拟合得到接近零频率的非线性水动力导数.该结论和方法可为供应船舶操纵性能的预报提供一定的参考.

[1]吴宝山.第七届船舶力学学术委员会操纵性学组工作总结报告[C].第七届船舶力学学术委员会全体会议论文集，2010:51-61.

[2]邹早建.船舶操纵性研究进展[C].第六届船舶力学学术委员会全体会议专集，2006:54-67.

[3]ITTC. The Manoeuvring Committee[C]. Proceedings of 25th ITTC-Volume I,2008:143-203.

[4]TURNOCK S R, PHILLIPS A B, Furlong M. Urans simulations of static drift and dynamic manouveres of the kvlcc2 tanker[C]. Proceedings of Workshop on Verification and Validation of Ship Maneuvering Simulation Methods (SIMMAN2008),2008.

[5]PAN Y C, ZHANG H X, ZHOU Q D. Numerical prediction of submarine hydrodynamic coefficients using CFD simulation[J]. Science Direct Journal of Hydrodynamics,2012,24(6):840-847.

[6]张赫，庞永杰，李晔.基于FLUENT软件模拟平面运动机构试验[J].系统仿真学报，2010，22(3):566-569.

[7]石爱国，闻虎，李理,等.船舶浅水水动力导数的数值计算[J].中国航海，2011，34(3):69-82.

[8]邹早建，杨勇，张晨曦.深浅水中KVLCC船体横荡运动水动力数值计算[J].水动力学研究与进展，2011，26(1):85-92.

[9]杨勇，邹早建，楼鹏宇.小振幅平面运动机构试验粘性流场数值模拟[C].第二十三届全国水动力学研讨会第十届全国水动力学学术会议文集，2011:190-195.

[10]杨勇.非定常操纵运动船体水动力数值计算[D].上海：上海交通大学，2011.

[11]邢磊.三体船水动力导数及操纵性能预报研究[D].哈尔滨：哈尔滨工程大学，2012.

[12]程捷，张志国，蒋奉兼,等.平面运动机构试验的数值模拟[J].水动力学研究与进展，2013，28(4):460-464.

[13]TOXOPEUS S L. Using CFD calculations to improve predictions of ship manoeuvres[J]. RINA Developments in Marine CFD,2011(3),13-25.

[14]NOBUAKI S, PABLO M, CARRICA, et al. URANS simulations of static and dynamic maneuvering for surface combatant: verification and validation for forces, moment, and hydrodynamic derivatives[J]. Journal of Marine Science and Technology,2012,17(4):422-445.

[15]李冬荔.粘性流场中船舶操纵水动力导数计算[J].哈尔滨工程大学学报，2010，31(4):421-427.

[16]刘山.基于CFD技术数值模拟平面运动机构试验[D].武汉：武汉理工大学，2012.

[17]林洪波.供应船的稳性特点[J].天津航海,1987(3):27-29.

[18]贾欣乐，杨盐生.船舶运动数学模型-机理建模与辨识建模[M].大连：大连海事大学出版社，1999.

[19]赵小仨，罗薇，柯枭冰，等.呆木对供应船线性水动力导数影响[J].船舶工程，2014，36(3):33-36.

[20]赵小仨，罗薇，柯枭冰，等.供应船纵向加速度导数的估算方法[J].中国航海，2014，37(4):88-91.

Numerical Simulation and Research of Nonlinear Hydrodynamic Derivatives

ZHAO Xiaosa1,2)XU Haixiang1,2)

(KeyLaboratoryofHighPerformanceShipTechnologyofMinistryofEducation,WuhanUniversityofTechnology,Wuhan430063,China)1)(Collegeoftransportation,WuhanUniversityofTechnology,Wuhan430063,China)2)

For the problem that the nonlinear dynamics of ships are difficult to define and forecast, this article uses CFD commercial software FLUENT in combination with dynamic mesh to simulate the large amplitude planar motion mechanism test for a supply vessel. The range of nonlinear motion amplitude of the supply vessel is obtained by comparing the pressure nephograms of different scenarios. Besides, the test scheme of the supply vessel’s nonlinear motion is designed. The large drift angle of oblique motion and large amplitude of pure surge, pure sway, pure yaw, combination movement with several frequencies are simulated, respectively. Finally, the nonlinear hydrodynamic derivatives close to zero frequency are calculated.

nonlinear hydrodynamic derivatives; large amplitude pmm test; numerical simulation; supply vessel

2016-12-14

*国家自然科学基金项目 (61301279, 51479158)、中央高校基本科研业务费专项资金项目(163102006)资助

U661.33

10.3963/j.issn.2095-3844.2017.01.014

赵小仨(1989—)：女，工学硕士，实验员，主要研究领域为船舶水动力研究