实船阻力及流场数值预报方法

易文彬，王永生，杨琼方，李剑

(海军工程大学 动力工程学院， 湖北 武汉 430033)

船阻力预报大部分是基于模型试验进行的，但是模型试验一般是在傅汝德数相等的条件下进行的，与实船雷诺数并不相等。在实际工程中，可以通过模型试验测得模型的总阻力，然后通过外推经验公式(二因次法或三因次法)和一些船体粗糙度等有关的经验修正项得到实船阻力。这些方法具有很强的工程实用性，但是包含了大量的经验修正项，船模与实船相关换算缺乏严格的理论支撑，也不能给出实船雷诺数下的流场信息。

目前基于CFD的船舶数值计算主要集中于模型尺度的雷诺数(107)，对实尺雷诺数(109)下的研究还比较少。但是随着船舶计算流体力学的发展和工程实践的需要，实尺度雷诺数下阻力和流场的精确计算获得了更加越来越多的关注[1]。Raven等[2]通过CFD的方法计算了模型和实船的阻力并且分析了船模实船阻力换算过程中各阻力成分的尺度效应。Bhushan等[3]通过Athena船模和实船的数值计算与校验，验证了一种新型的壁面函数模型。其中实尺计算值考虑了粗糙度的影响，与船模试验外推值吻合良好。刘志华等[4]提出了一种基于设置流体运动粘度系数实现雷诺数相似的计算理论，能够较为快速准确地求取实船雷诺数下的自航因子。倪崇本等[5]提出将实尺度下势流理论与模型尺度下湍流理论相结合的方法来求取实船的阻力，虽然阻力预报有较高的精度，但是不能给出实船雷诺数下的流场信息。

本文采用基于数值计算的3种方法预报了实船阻力、波形及流场信息，对计算结果进行了比较和分析。在虚流体粘度方法中分析了粗糙度对各阻力成分的影响，比较了阻力与流场的计算结果与其他方法计算结果的差异。

1 基于CFD的船舶总阻力预报方法

1.1 数值计算模型及控制方程

采用有限体积法离散控制方程，不可压缩粘性流体的控制方程如下：

(1)

(2)

式中：ui分别为流体质点在i方向的速度分量，fi是质量力，p为流体的压力，μ是相体积分数平均的动力粘度系数，μt为湍流动力粘性系数。

上述方程需要结合湍流模型封闭方程组，常见的两方程湍流模型有k-ε、RNGk-ε、k-ω、k-ωSST、RSM等。本文选取的是k-ωSST湍流模型。

1.2 船舶阻力试验换算方法

船模阻力试验换算方法最早由傅汝德提出，该方法将船体总阻力分成摩擦阻力Rf和剩余阻力Rr，并认为摩擦阻力仅与雷诺数Re有关，剩余阻力仅与傅汝德数Fr有关，且两者互不干扰[6]。因此总阻力系数可以表示为

Ct(Fr,Re)=Cf(Re)+Cr(Fr)

(3)

因此实船尺度下的阻力系数可由船模阻力试验中测得的阻力系数表示为

(4)

式中：Cf0代表相当平板摩擦阻力系数，Cr为剩余阻力系数，下标s和m分别代表变量对应于实船和模型。在此基础上，休斯提出了三因次换算方法，引入形状因子1+k来表示船体粘性阻力和相当平板摩擦阻力之间的关系：

(5)

(6)

引入形状因子的三因次法较二因次法合理，其中形状因子可由低速拖曳或叠模试验确定。考虑到实船船体表面粗糙度的作用以及船模实船阻力换算过程中由于雷诺数不相等导致的尺度效应，需要补贴一定的经验系数，即船模实船换算补贴系数Ca，Ca可由下式确定：

(7)

式中：Lpp为垂线间长，ks为粗糙度表观高度，一般可取ks=0.15 mm。

1.3 船模实船雷诺数和傅汝德数全相似条件

为保证船模和实船总阻力系数相等，必须满足实船和船模的雷诺数Re和傅汝德数Fr相等：

(8)

(9)

式中：L为船长，υ为流体运动粘度，V为来流速度，g为重力加速度。在船模试验时不能满足全相似定律，因此工程实践中船模试验都是在傅汝德数相等的条件下进行的。为实现实船和船模的雷诺数和傅汝德数全相似，在数值计算中可令

(10)

此时，式(8)、(9)同时成立。引入虚拟流体运动粘度系数之后，流体仍为牛顿流体，流场的基本运动规律不变，因此湍流模型及控制方程不需要调整。

2 计算结果及分析

2.1 基于模型尺度数值模拟外推方法

本文的主要研究对象是DTMB 5415，该船型是ITTC推荐的军舰类型中唯一有大量公开试验数据的平台。美国的DTMB，IIHR以及意大利的INSEAN三家研究机构对该船型进行了全面的试验。本文对DTMB 5415进行了流场模拟和阻力计算，并以试验数据[7-8]为校验。因DTMB 5415没有相关的实船数据，本文暂将模型放大20倍作为实船。

计算区域入口取船艏向前延伸1倍船长处，出口取船艉向后延伸2倍船长处，侧边界及下方边界均取1倍船长，上方边界取水线上方0.5倍船长。采用VOF的方法来追踪自由液面，湍流模型采用k-ωSST模型，对流项离散选用二阶离散格式。本文采用分块结构化网格，船体周围采用O型网格以便捕捉船体附近的边界层，同时对船艏、船艉及静水面进行网格加密，网格节点数320万。满足计算要求。

Fr为0.15、0.28、0.41 共3个航速下阻力计算值与试验值如表1所示。其中Fr=0.28及Fr=0.41航速下考虑了船体姿态的变化。从表1可以看出：通过数值模拟能够较好地预报船模的阻力，阻力计算值误差在3%以内。若采用三因次换算方法求取实船阻力，还必须求得形状因子1+k。ITTC建议用傅汝德数低于0.15的船模阻力试验数据来求取形状因子。本文通过叠模计算求取傅汝德数在0.1～0.15时的形状因子，计算结果如表2。

表1 DTMB 5415船模总阻力系数

表2 DTMB 5415叠模阻力计算

取多个速度下形状因子的平均值，1+k=1.083。由低速时拖曳船模试验数据通过普鲁哈斯卡方法换算得到1+k=1.15。叠模计算得到的形状因子1+k偏小，因为船模低速拖曳试验时兴波虽然不明显，但是仍然存在兴波阻力成分，所以得到的形状因子比叠模计算值大。考虑船模实船换算补贴Ca(按照式(7)本船取为0.000 509)，将模型试验与数值计算得到的船模阻力系数Ctm按照三因次法(1+k分别取1.15和1.083)外推得到的实船阻力系数Cts。其值如表3所示。

表3 1+k法预报的阻力

其中相对偏差定义为

(11)

基于模型尺度数值模拟外推方法得到的实船阻力与由模型试验外推得到的阻力差别在7%以内。因叠模计算得到的形状因子1+k与通过低速拖曳试验得到的1+k相比偏小，如果采用相同的船模实船换算补贴，通过计算预报的实船阻力与通过试验预报的实船阻力存在较大的偏差。

2.2 “虚流体粘度”方法

为实现船模、实船雷诺数Re和傅汝德数Fr全相似条件，在对船模进行数值模拟时将水的运动粘度按照全相似的条件设置为一个虚拟的粘度系数(本文将水的运动粘度系数取为1.067×10-8m2/s)。文献[9]表明，在船模雷诺数下，如果要准确模拟船模的摩擦阻力，船体表面第一层网格的厚度要满足y+在30～ 300之间，但是在实船雷诺数下，y+最大值可能远远大于1 000。于是本文基于船体摩擦阻力变化对船体表面第一层网格的厚度做了相关的研究：本文针对Fr=0.28航速下，划分了5套网格，船体表面第一层网格厚度分别为0.2、1、2、3、5 mm，计算得到的摩擦阻力系数如表4。

表4 不同网格下摩擦阻力计算结果

从表4可以看出，针对本算例情形，船体表面第一层网格厚度应该使y+在5 000～ 10 000时能够准确模拟实船雷诺数下的摩擦阻力。进行实船阻力计算，还必须要考虑粗糙度的影响。目前在数值计算中，一般是采用平均沙粒粗糙度来代替船舶实际粗糙度的影响，文献[10-11]表明将实船平均沙粒粗糙度取为hs=0.03 mm与粗糙度表观高度ks=0.15 mm是等效的。若船模实船相对粗糙度相等，则

hs/Ls=hm/Lm

(12)

将船模平均沙粒粗糙度hm取为0.001 5 mm。光船阻力系数和考虑粗糙度计算得到的实船阻力系数如表5。定义相对偏差为

(13)

表5 虚流体粘度方法阻力预报

从表5可以看出：考虑粗糙度的虚流体粘度系数方法阻力计算值与船模试验外推值在傅汝德数从0.15～0.41航速范围内相差不到4.5%，表明这种方法能够较好的预报实船的阻力。在计算中同时考虑粗糙度和实船雷诺数的作用，克服了船模实船换算的尺度效应，可以直接求得实船的阻力系数。此方法计算所需的网格与普通模型数值模拟的网格量接近，且不需要进行叠模计算求取形状因子，可以快速准确地实现实船阻力预报。

为进一步分析粗糙度的影响，将考虑粗糙度前后各阻力成分的变化用图1，图2表示如下。

图1 粗糙度对摩擦阻力的影响

图2 粗糙度对压阻力影响

不考虑粗糙度时，由虚流体粘度方法计算得到的摩擦阻力系数与相当平板摩擦阻力系数是很接近的，压阻力系数与模型雷诺数下计算得到的压阻力系数差别不大。考虑粗糙度之后，发现摩擦阻力系数增大了12%～31%，压阻力系数基本不变，表明粗糙度对摩擦阻力影响较明显，对压阻力影响不大。

2.3 实尺度船舶RANS计算

针对实尺度计算，采用与模型尺度下相似的控制域，并且对船体和水线面附近的网格进行加密，网格数目为1 500万。为缩短计算时间，本文以层流的方法得到的结果作为高精度湍流计算的初值，以光滑船体的计算结果作为粗糙船体计算的初值，实践证明，按照这种方法可以极大地缩短计算时间。实船的平均沙粒粗糙度取0.03 mm。不考虑粗糙度和考虑粗糙度的实船阻力计算结果如表6。

表6 实尺度计算阻力预报

考虑粗糙度的实船阻力CFD计算值与船模试验外推值偏差在3.5%以内，表明实船的数值模拟能够提供较精确的实船阻力。实尺度计算不需要相关的经验系数补贴，能够克服船模实船换算的尺度效应，提供实尺度下的流场信息。随着计算机技术的发展，实尺度下的数值模拟表现了越来越强的工程应用意义。

2.4 波形和流场比较和分析

图3为Fr=0.28航速下船体的兴波，船模的数值模拟与船模试验值的波形基本吻合。数值计算预报的波形比较接近，但是波幅有细微的差别，其中实船雷诺数下船艉的波幅更大。

图3 船体附近兴波

波高的比较结果见图4。定义原点为首垂线与水线的交点，x轴指向船艉，y轴指向船体右侧，z轴指向船体上方。

图4 Y/Lpp=0.172处的波高

通过y/Lpp=0.172处波高的对比，可以看出：

1)模型试验与船模数值计算得到波峰和波谷的位置吻合较好，但是计算得到的首波峰较小，可能是数值耗散的原因。在1.5Lpp以后由于网格较疏，波高与试验值有一定的偏离。

2)虚流体粘度给出的结果与实尺度计算得到的波高吻合较好，表明虚流体粘度的方法可以较好地预报实船的波高。

3)波高在0～Lpp处都吻合较好，但是在船艉处实船雷诺数下波幅要大于模型雷诺数下的波幅，这与之前波形的分析结果也是一致的。因为在模型雷诺数下，船艉处的边界层厚度相对较大，粘性对兴波的影响更加明显，导致波幅较小。

图5为桨盘面处伴流场，可以看出：

1)模型尺度数值模拟能够和试验数据较好的吻合，但是在一些数据点处仍有偏离。

2)虚粘度方法与实船计算给出了较为一致的结果，且伴流作用较模型雷诺数下弱。主要原因是模型雷诺数下船尾的边界层较厚，对伴流的影响更加显著。

图5 桨盘处纵向伴流(X/Lpp=0.935,Z/Lpp=-0.02)

通过对波形、波高和伴流场的比较可以看出：虚流体粘度方法和实尺度计算得出的结果较为一致且体现出了实船雷诺数下流场特征；粘性的作用相对较弱，导致伴流较弱，兴波更强。

3 结论

本文基于CFD对DTMB 5415实船尺度下的阻力和流场进行了计算和预报。通过对计算结果的进一步分析，得到以下结论：

1)通过合理地设置船体表面第一层网格厚度和粗糙度，可以通过虚流体粘度方法或实尺度计算准确地求取实船的阻力。

2)虚流体粘度方法和实尺度计算均能直接求取实尺度雷诺数下的阻力和流场，但虚流体粘度方法网格量少，计算量较小，表现了较强的工程意义。

3)基于船模尺度数值模拟外推值与船模试验外推值偏差较大，主要原因是形状因子的计算值与试验采用的数值差别较大。有关形状因子的计算有待进一步探讨。

随着尺度效应和粗糙度相关研究和试验的开展，本文的研究工作有待进一步深入。

参考文献：

[1]STERN F, YANG Jianming, WANG Zhaoyuan, et al. Com-

putational ship hydrodynamics: nowadays and way forward[C]//29th Symposium on Nnaval Hhydrodynamics. Gothenburg,Sweden, 2012:26-31.

[2]RAVEN H C, VAN D, STARKE A R, et al. Towards a CFD-based prediction of ship performance—Progress in predicting full-scale resistance and scale effects[J]. International Journal of Maritime Engineering, 2009,150(5):14-18.

[3]BHUSHAN S, XING T, CARRICA P, et al. Model- and full-scale URANS simulation of Athena resistance, powering, seakeeping, and 5415 maneuvering[J]. Journal of Ship Research,2009,53(4): 179-198.

[4]刘志华, 熊鹰, 韩宝玉. 雷诺相似船模预报实船推进因子的数值方法[J]. 哈尔滨工程大学学报,2008, 29(7):658-662.

LIU Zhihua, XIONG Ying, HAN Baoyu. Accurately predicting the propulsion factors of a ship using a Reynolds-similarity ship model[J]. Journal of Harbin Engineering University, 2008,29(7):658-662.

[5]倪崇本, 朱仁传, 繆国平, 等.一种基于CFD的船舶总阻力预报方法[J].水动力学研究与进展,A辑, 2010, 25(5):579-586.

NI Chongben, ZHU Renchuan, MIAO Guoping, et al. A method for ship resistance prediction based on CFD computation[J]. Chinese Journal of Hydrodynamics, 2010,25(5):579-586.

[6]盛振邦, 刘应中. 船舶原理[M]. 上海：上海交通大学出版社,2003:154-157.

[7]OLIVIERI A, PISTANI F, AVANZINI A, et al. Towing tank experiments of resistance, sinkage and trim, boundary layer, wake, and free surface flow around a naval combatant INSEAN 2340 model[R].Iowa:The University of Iowa,2005.

[8]STERN F, LONGO J, PENNA R, et al. International collaboration on benchmark CFD validation data for surface combatant DTMB model 5415[C]//Proceedings of 23rd ONR symposium on naval hydrodynamics. Iowa,2000:17-22.

[9]ANSYS Inc. ANSYS CFX-solver theory guide[M]. (s.l.):ANSYS Inc, 2007.

[10]HE Chengbin, HUANG Jianguo, HAN Jin, et al. Cyclic shift keying spread spectrum underwater acoustic communication[J]. Acta Phys Sin, 2009, 58(12): 8379-8385.

[11]SCHULTZ M P. Effects of coating roughness and biofouling on ship resistance and powering[J]. Biofouling, 2007,23(13):331-341.

[12]ECA L, HOEKSTRA M, RAVEN H C. Quantifying roughness effects by ship viscous flow calculations[C]// 28th Symposium on Naval Hydrodynamics. California, 2010:12-17.