基于CFD的船舶横摇数值模拟与粘性效应分析

罗天，万德成

1上海交通大学船舶海洋与建筑工程学院，上海200240

2上海交通大学海洋工程国家重点实验室，上海200240

3高新船舶与深海开发装备协同创新中心，上海200240

基于CFD的船舶横摇数值模拟与粘性效应分析

罗天1，2，3，万德成1，2，3

1上海交通大学船舶海洋与建筑工程学院，上海200240

2上海交通大学海洋工程国家重点实验室，上海200240

3高新船舶与深海开发装备协同创新中心，上海200240

［目的］船舶横摇的准确预报对于船舶耐波性、稳性以及操纵性的研究具有十分重要的意义，但船舶横摇运动受流场粘性效应的影响很大，计算中存在较多非线性因素，因而并不适于常用于研究船舶运动的传统势流理论。为解决这一问题，［方法］采用基于OpenFOAM软件开发的naoe-FOAM-SJTU求解器，分别通过欧拉方法和RANS方法对S60船模典型二维横剖面的强迫横摇进行数值模拟，同时，模拟分析DTMB 5512不同三维船模的自由与强迫横摇运动。［结果］在成功将横摇阻尼力矩的不同成分分别计算出来后发现，其中漩涡阻尼所占的比例最大，摩擦阻尼所占的比例最小，而舭龙骨在一定的横摇角度范围内则减小了横摇阻尼力矩。［结论］该结果揭示了横摇参数对船舶横摇运动及横摇阻尼力矩的影响，对准确预报船舶的横摇运动具有重要意义。

船舶；横摇运动；阻尼系数；阻尼力矩；RANS方法；欧拉方法

0 引 言

在现代船舶设计中，船舶横摇运动的准确预报十分重要，然而流场的粘性效应对船舶横摇运动影响较大，使得船舶横摇运动计算中含有很多非线性因素，因此传统的势流理论无法准确预报船舶横摇运动。

船舶横摇阻尼系数是决定船舶横摇预报准确性的关键因素。船舶横摇阻尼系数通常可以分为5个部分［1］：摩擦阻尼系数、波浪阻尼系数、漩涡阻尼系数、升力阻尼系数和附体阻尼系数。这种分类便于逐个研究外界因素对阻尼系数不同成分的影响，但在该分类下，5种阻尼系数之间是假定线性相关的。对于这5个部分，Ikeda等［2-3］于1977～1978年根据船模试验提出了相对应的经验公式，但船模试验的效率不高，经验公式的适用范围有限，因此这2种方法并不能大规模地用于船模横摇研究。

如今，随着计算机技术的飞速发展，计算流体力学［4］（Computational Fluid Dynamics，CFD）被广泛应用于船舶相关研究。在船舶横摇运动模拟分析方面，学者们相继进行了研究。Luca与Stefano［5］利用CFD软件针对相似横剖面的单体船、双体船和小水线面船进行了横摇运动模拟；Miller等［6］应用RANS（Reynolds Averaged Navier Stokes）方法对存在舭龙骨的三维圆柱体进行了预报；Wilson等［7］应用CFDSHIP-IOWA求解器预报了船舶在不同船速下的自由横摇以及强迫横摇；Yang等［8］针对S60船模的自由与强迫横摇运动进行了模拟并计算了其阻尼系数。但以上研究都只得到了总的阻尼系数，或是仅将摩擦阻尼系数分离出来进行了分析。

本文将应用课题组自行开发的naoe-FOAMSJTU求解器［9］模拟船模运动与流场的变化。该求解器是基于开源软件OpenFOAM［10］（Open Field Operation and Manipulation）而开发，以RANS方程为控制方程，采用有限体积法（Finite Volume Method，FVM）离散控制方程，流体体积法（Volume of Fluid，VOF）技术捕捉自由液面，用SSTk-ω湍流模型封闭方程，应用压力的隐式算子分割算法（Pressure-Implicit Split-Operator，PISO）［11］耦合求解速度和压力。而对于船体的运动，则采用动态网格技术以及六自由度运动求解器来模拟。naoe-FOAM-SJTU求解器已应用于船舶工程、近岸工程、海洋工程等各类复杂流动问题的数值模拟，包括各类水波的数值造波与消波、船桨舵匹配自航推进与操纵运动、节能装置水动力性能、船舶运动与粘性阻尼、上浪砰击、浮式平台运动与气隙、液化天然气（LNG）液舱晃荡、Spar平台涡激运动等复杂问题，相关计算结果也得到了很多试验的证实。运用该求解器，吴建威等［12］成功预报了螺旋桨敞水性能，Shen等［13］得到了一种非规则波的生成方法，李鹏飞等［14］基于致动线模型数值模拟了风力机尾流场，尹崇宏等［15］模拟了实尺度VLCC船的阻力问题，Cao等［16］数值模拟了多方向海洋波问题，Zha等［17］数值模拟了双体船兴波阻力问题。Shen等［18］还数值模拟了船舶自航运动和船—桨—舵的相互作用问题。以上数值模拟结果均表明naoe-FOAM-SJTU求解器在解决船舶与海洋工程相关问题方面是可靠的。

本文将分别采用欧拉方法与RANS方法模拟二维船模横剖面的强迫横摇运动，并从总阻尼系数中分离出波浪阻尼系数、摩擦阻尼系数以及漩涡阻尼系数。对于三维船模，拟通过比较带舭龙骨的船体与裸船体的横摇运动，从总阻尼系数中分离出摩擦阻尼系数与舭龙骨阻尼系数。由于本文中所涉及到的横摇运动并未考虑船模横荡运动的影响，因此升力阻尼系数不在分析范围内。

1 数值方法

1.1 RANS方法

对于naoe-FOAM-SJTU求解器，以RANS方程为控制方程：

式中：x为网格节点的位置；t为时间；U为速度场；Ug为网格移动速度；pd=p-ρg·x，为动压力；ρ为液体或者气体的密度；g为重力加速度向量；fσ为表面张力项；fs为源项；μeff=ρ(ν-νt)，为有效动力粘度，其中ν为运动粘度，νt为涡黏度。

湍流模型SSTk-ω用来封闭该控制方程，边界可压缩的VOF方法用来追踪和捕捉自由液面的变化，其输运方程定义如下：

式中：Ur为用于压缩界面的速度场；α为积分数，代表网格单元中流体所占的体积比，其值介于0和1之间：

控制方程采用FVM法进行离散，即将计算域离散成一系列小单元，计算流场信息存储在网格单元中心，最终根据单元中心的值插值得到单元面的值。对于离散后得到的压力速度耦合方程，则采用PISO算法进行迭代求解，其主要由分预测—修正—再修正三步组成。具体步骤与算法详见文献［9］。

1.2 欧拉方法

本文采用欧拉方法对二维S60船模横剖面的强迫横摇进行计算。由于摩擦阻尼系数和漩涡阻尼系数均与流场的粘性相关，同时，本文假定波浪阻尼系数与流场粘性无关，因此由欧拉方法最终得到的为波浪阻尼系数。

欧拉方法的实现同样是基于naoe-FOAMSJTU求解器。可将求解器中粘性项的数值设置成0，并改变流场的边界条件来实现欧拉方法，具体同第1.1节。

1.3 静水阻力系数

为了验证欧拉方法的可行性，同时验证波浪与粘性无关，本文运用欧拉方法以及RANS方法对KCS船模在静水中的运动进行模拟，并求解静水阻力系数。

采用RANS方法可以直接计算得到粘性流场中船舶的总阻力，从而求得阻力系数

式中：F为船舶总阻力；S为船舶湿表面积；V为船速。

而采用欧拉方法则仅能计算出兴波阻力系数Cw，摩擦阻尼系数Cf采用1957ITTC摩擦公式：

式中，雷诺数Re=VL/ν，其中L为船长。最终的阻力系数C=Cf+Cw。

1.4 横摇阻尼系数

本文中，为了与试验数据进行对比验证，需要计算相关的横摇阻尼系数。强迫横摇与自由横摇阻尼系数的计算方法并不一致，具体如下。

1.4.1 强迫横摇

对于强迫横摇，横摇角度随时间的变化可以

表示为如下方程：

式中：θ0为初始横摇角；ω为横摇频率。

横摇阻尼力矩可以表示为

式中，Beq为阻尼系数。假设回复力系数C(θ，t)是关于θ的多项式，则

当θ=0时，可以推导得到

式中：M0为船舶横摇角度为0°时船舶所受到的阻尼力矩；B为船宽。

由以上方程可以看出，对于某一种确定的工况，阻尼系数与横摇阻尼力矩相关，因此，横摇阻尼力矩相应地也可以分为5个部分。为直观分析，本文将优先从横摇阻尼力矩曲线以及最大力矩来进行分析。

本文采用naoe-FOAM-SJTU求解器求得的总阻尼力矩可以分为剪切和法向2个部分。剪切部分是由船体周围的剪切力所产生，即摩擦力产生的阻尼力矩，相对应的为摩擦阻尼系数。而法向部分则是由船体周围的不平衡压力所造成，对应于波浪阻尼系数和漩涡阻尼系数，其中假定波浪阻尼系数与粘性无关，是采用欧拉方法计算得到。舭龙骨的存在增加了船体摩擦，同时也改变了船体周围的涡量分布和涡量大小，因此对于剪切与法向2个部分均有影响。4个部分（漩涡、波浪、摩擦、附体）的横摇阻尼力矩的关系如下：

式中：Mp和Mf为采用RANS方法求得的阻尼力矩法向部分与剪切部分，其中Mp是由漩涡、波浪与舭龙骨共同影响产生，Mf是由摩擦部分产生；Mpw为采用欧拉方法求得的阻尼力矩法向部分，是由波浪影响产生的；Mpe为由漩涡影响产生的横摇力矩；Mpb为由舭龙骨单独影响产生的横摇力矩的法向部分；Mff为由裸船体的摩擦产生的横摇剪切力矩；Mfb为由舭龙骨单独的摩擦产生的横摇剪切力矩。

对于自由横摇运动，除了需要通过横摇角度的时历曲线来估算横摇衰减周期外，还需通过该曲线来计算船舶自由横摇的阻尼系数 μϕϕ［19］。图1所示为典型的船舶自由衰减横摇曲线。

根据船舶耐波性理论［20］，由该曲线可以得到相邻两振幅的差Δφ=φk-φk+1以及均值φm=(φk+φk+1)/2，以 φm为横轴、Δφ为纵轴作消灭曲线如图2所示。

则船舶自由横摇阻尼系数可通过下式求得：

2 计算模型和网格

本文一共使用了3种船型：KCS船型、S60船型以及DTMB 5512船型。计算网格由OpenFOAM自带的网格生成工具SnappyHexMesh生成。为了提高计算精度，捕捉自由面兴波以及船体附近的涡量变化情况，对自由面和船体壁面附近处的网格进行了适当加密。除此外，为了实现船模的横摇运动，还在计算中引入了动态网格（Dynamic mesh）技术。

2.1 KCS船型

在通过静水阻力验证欧拉方法的可行性时，为便于与试验数据［21］进行对比，采用了与试验相同的船型——KCS船型，其船模参数如表1所示。

KCS船的船体视图如图3所示。

由于船体的对称性，仅在船模中纵剖面的右舷侧划分了网格，计算域尺寸为：-1.0L≤x≤4.0L，0≤y≤1.5L，-1.0L≤z≤1.0L。网格总量为108×104，船体艏、艉网格以及背景网格如图4所示。

（1）通过对国内外现有软件的应用分析调研，确定基于BIM技术的绿色建筑规划设计软件架构和数据结构，建立基于BIM技术的三维建筑模型

2.2 S60船型

对于本文中讨论的典型的二维船型横剖面，为便于与试验数据［22］进行对比，采用了与试验相同的船型及剖面——S60船型的典型横剖面，船型主要参数如表2所示。

对于该二维横剖面的横摇模拟，横摇轴为x轴，横摇的轴心O为船舶横剖面中心线与船模横剖面静止时水线的交点。计算域尺寸为：-2.5B≤y≤2.5B，-2B≤z≤2B。另外，还对该横摇算例进行了网格收敛性分析，网格量分别为4×104，7×104和14×104，计算网格如图5所示。

2.3 DTMB 5512船型

对于本文中讨论的典型三维船型，为便于与试验数据进行对比，采用了与试验相同的船型——DTMB 5512船型，船型主要参数如表3所示。

该船体的裸船模型以及带舭龙骨的船模如图6所示。

对于横摇运动，该船型的横摇中心与重心重合，横摇轴为中纵剖面上通过重心且平行于水线面的直线。计算域尺寸为：-1.0L≤x≤3.0L，-2.0L≤y≤2.0L，-1.0L≤z≤0.5L。为了验证网格的可靠性，该算例也做了网格收敛性验证，网格量分别为97.8×104，232×104和646×104，其中，232×104的计算网格如图7所示。

3 验证与计算结果

3.1 S60典型二维横剖面横摇计算

本文结合欧拉以及RANS这2种方法对S60二维横剖面的强迫横摇运动进行了模拟计算，最终将横摇阻尼力矩分为了波浪、摩擦和漩涡3个部分。

3.1.1 欧拉方法可行性验证

为了验证欧拉方法的可行性以及其对船体周围波浪预报的准确性，本文使用naoe-FOAMSJTU求解器，采用RANS和欧拉2种方法对不同傅汝德数Fr下的船舶静水阻力进行了计算并与试验数据进行了对比，计算结果如表4所示。

从表4的结果可以看出，由该求解器应用欧拉方法得到的结果与试验结果相近，即欧拉方法计算得到的兴波阻尼部分是可信的，可以由该方法计算横摇阻尼力矩中的波浪阻尼力矩。

当Fr=0.27时，由2种方法计算得到的波形比较如图8所示。图中，上半部分为欧拉计算结果，下半部分为RANS计算结果。

由上述阻力系数比较可以看出，采用欧拉方法和RANS方法所得结果与试验结果的误差在允许的范围之内，说明欧拉方法对于非粘性流场的计算结果是可靠的。从采用欧拉方法与RANS方法模拟得到的稳定波形图上来看，两者在船体附近的波形图基本相同，因此可以推断流场的粘性效应对于船体附近的兴波并没有太大的影响，继而可以推论，在横摇模拟中，采用欧拉方法计算波浪阻尼部分的做法是可行的。

3.1.2 RANS方法可行性验证

为了与Ikeda等［2-3］的试验数据进行对比，本文采用RANS方法在相同的横摇工况下对S60典型二维横剖面进行了强迫横摇计算，横摇工况如表5所示。

首先选用0.17 rad进行收敛性验证，结果如表6所示。

由表6中数据可知，网格2与网格3的误差均在允许范围内，但是考虑到效率问题，在后面的计算中，将选择网格2进行计算。

为了进一步验证网格2计算结果的可靠性，使用网格2对另外3种横摇角度下的强迫横摇进行了计算模拟，计算结果如表7所示。

由表7中的数值结果可以看出，误差均在允许范围内，可以证实采用RANS方法以及网格2计算该二维船模横剖面强迫横摇运动所得结果是可以信赖的。

3.1.3 S60典型二维横剖面强迫横摇计算结果

为了单独分析二维船模横剖面强迫运动时横摇阻尼系数的不同成分，本文运用欧拉方法以及RANS方法对于同一种横摇工况下的船舶横摇进行了模拟，具体如表8所示。

在该工况下，二维船模横剖面强迫横摇的力矩时历曲线如图9所示。

图9中，黑色曲线为采用欧拉方法算得的船舶所受力矩的垂向部分 Mpw，红色曲线为采用RANS方法算得的船舶所受摩擦力矩Mff，蓝色曲线为采用RANS方法求得的力矩垂向部分与Mpw之差，即Mpe。

粘性流场中，船舶强迫横摇时船体横剖面周围流场的涡量分布如图10所示。

图10为第（n-1）个周期内4张典型时刻的涡量分布图，其中T表示强迫横摇周期，t表示该图所对应的时刻。船舶强迫横摇时粘性流场与非粘性流场中的涡量分布图如图11所示。

由图9可以看出，在强迫横摇阻尼力矩中，漩涡影响部分相较于其他部分占比要大得多，因为在船体横摇过程中周围会产生很多漩涡，造成船体周围压力的不均衡，相比于摩擦以及由波浪引起的力矩，这部分产生的力矩影响更大。同时，这种结果与Ikeda等［2-3］的试验结果和Korpus等［23］的数值计算结果均很相似。另外值得注意的是，波浪部分产生的力矩与其他2个部分力矩的方向是相反的，其原因分析如图12所示。

由图12可知，波浪对船舶形成的力矩使船趋向于倾覆方向，因此与另外2个部分力矩的方向是相反的，但因其所占总力矩的比例较小，故总的横摇力矩还是与漩涡产生的横摇力矩方向相同。

另外，从力矩时历曲线图（图9）中可以看出，漩涡阻尼部分和摩擦阻尼部分达到峰值与零点是有时间差的，摩擦阻尼部分出现峰值的时间比漩涡阻尼部分要迟。船模所受到的漩涡阻尼部分由二维船模横剖面横摇时流场的具体状态所决定，由上面的数据可以看出，峰值出现在横摇角度最大的时刻，对应于计算结果的涡量图，发现此时船模横剖面周围的漩涡最多，尤其是舭部，造成了船舶周围压力的极度不平衡，因此，力矩的漩涡阻尼部分达到峰值。而摩擦阻尼部分是由船模横剖面临近的液体与船模表面的相对速度所决定，当船模横剖面处于横摇角度最大时，船模表面的速度为0，而临近的流体由于惯性仍有速度，之后，船模向反方向运动，使得相对速度继续增大，即力矩的摩擦阻尼部分会继续增大，因此摩擦阻尼部分峰值出现的时间会比漩涡阻尼部分晚。

分析船舶横摇涡量图（图10）可知，漩涡产生于二维船模横剖面的壁上，在船模横摇角度从最大开始逐渐变小时，由于相对速度在船模周围开始形成一层很薄的、方向一定的漩涡，因此随着船舶逐渐左倾（或右倾），在船模横摇角度到另一个方向达到最大时，该方向上的漩涡已全部产生并附在船体表面形成较厚的一层。之后，船模反方向转动，立刻会在船模周围形成反方向的漩涡，中和一部分外层不同方向的漩涡，剩下没有被中和的部分则沿着船体底部和船侧方向分别排向船底与自由液面，在船模周围形成4处漩涡中心。这4处漩涡的方向由最初船体转动的方向所决定，当船体初始转动方向变动时，漩涡的方向将变为反向。

根据欧拉方法和RANS方法计算所得的无粘性与粘性流场中的压力与涡量分布图（图11）可以明显看出，当无粘性时，流场中并未产生明显的涡量，因此可以证实阻尼系数的漩涡部分仅由粘性影响产生，即漩涡阻尼系数是根据法向力矩的粘性部分而求得。

3.2 DTMB 5512三维船模横摇计算

本文采用RANS方法对DTMB 5512三维船型的横摇运动进行了模拟计算，最终将摩擦阻尼部分以及舭龙骨阻尼部分从总阻尼中分离了出来。

3.2.1 三维船模自由横摇计算验证

对于DTMB 5512三维船模，仅有在稳定航速下的自由横摇试验结果，因此本节采用同样的计算工况来模拟该船模的自由横摇运动以验证RANS方法和网格的可行性。计算工况具体如表9所示。收敛性验证结果如表10所示。

从表10中的数据可知，网格2与网格3的误差均在允许范围内，但是考虑到效率问题，在后面的计算中将选择网格2进行计算。

为了进一步验证使用该方法以及该网格所得结果的可靠性，使用网格2，根据试验数据采用了2种来流速度进行计算验证，2种工况具体如表11所示。

计算结果如图13和表12所示。由上面的结果可以看出，本文计算与试验间的误差均在6%以内，可以证明采用RANS方法以及网格2计算该三维船模横剖面强迫横摇运动结果可信。此外，从图13中可以发现，数值模拟结果与试验结果横摇周期的相对误差一定，但绝对误差值会随时间累积，在自由横摇的横摇角度时历曲线上，就表现为数值模拟结果与试验结果间的误差是随时间逐步增大的。

3.2.2 DTMB 5512三维船模强迫横摇计算结果

为了单独分析三维船模强迫运动时横摇阻尼系数的不同成分，采用RANS方法对同一工况下带舭龙骨的船模与裸船模进行了模拟计算，工况具体如表13所示。

在该工况下，强迫横摇的力矩时历曲线如图14所示。

图14中，黑色线为带舭龙骨船模所受的总力矩M，绿色曲线为裸船模所受力矩总力矩的垂向部分 Mpw+Mpe，蓝色曲线为裸船模所受的摩擦力矩 Mff，红色曲线为受舭龙骨影响的力矩的垂向部分Mpb，粉色曲线为受舭龙骨影响的力矩的切向部分Mfb。

粘性流场中船舶强迫横摇时船体横剖面周围流场的涡量分布如图15所示。

由图15中结果可以看出，在三维船模的强迫横摇中，与二维类似，仍然是Mpw+Mpe占总横摇阻尼力矩的大部分，而裸船模的摩擦力矩和带舭龙骨船模的摩擦力矩所占的比例则非常小。此外值得注意的是，对于带舭龙骨的船模，计算得到的横摇总阻尼力矩减小了，其具体原因如图16所示。

由舭龙骨造成的阻尼力矩变化的原因来源于2个方面：一是舭龙骨增加了船模表面积从而增加了摩擦阻尼；二是舭龙骨在船舶强迫横摇中对流场中涡量的数值以及分布有一定的影响从而改变了船舶受到的横摇阻尼力矩。其中，摩擦阻尼占据很小的一部分，在此阶段的分析中可以忽略不计。第2个因素首先考虑的是由舭龙骨附近产生的压力不平衡所造成的阻尼力矩。如图16所示，在船模转到最高点往平衡位置转的时候，总阻尼力矩是有利于船模转回平衡位置的，但此时由舭龙骨受到的流场的力造成的横摇阻尼力矩则是阻止船模转回平衡位置的。其次，由图15所示的涡量分布图可以看出，在有舭龙骨的情况下，在舭龙骨的范围内船体舭部形成了2个漩涡中心，改变了船模周围流场涡量的数值与分布。在这2种因素的综合影响下，流场对船模造成的 Mp在一定的横摇角度范围内减小了，其中也包括Mp的最大值。

与二维船模横剖面强迫横摇类似，三维船模的摩擦阻尼力矩时历曲线与阻尼力矩垂向部分的时历曲线具有一定的相位差。Mp与船体周围的涡量有关，在船模横摇达到最大角度时，船模附近的涡量分布最不均匀，涡量数值最大，因此Mp也达到最大值。当船模从横剖面处于横摇角度最大处向反方向运动时，临近流体与船模表面的相对速度会继续增大，即力矩的剪切力部分会继续增大，因此两者的时历曲线会有一定的相位差。

由图15所示的涡量图还可看出，船模附近漩涡以及漩涡中心形成的原理也与二维船模横剖面类似。漩涡在船模横摇时从船模表面开始形成，随着船模的运动，部分与新形成的漩涡相互中和消散，部分受船模运动的影响而向外移动。由于算例中的强迫横摇角度偏小，裸船模最终仅形成了2个漩涡中心，而带舭龙骨的船模在舭龙骨附近也会形成涡中心。图15还展示了船模横摇时自由表面的变化，可以看出，带舭龙骨船模与裸船模形成的行波基本相似。

4 结 语

本文采用RANS方法对S60船模二维横剖面的横摇进行模拟计算，验证了网格以及计算方法的可行性。在此基础上，结合欧拉方法下同一模型与工况的计算结果，将船舶阻尼系数分为了摩擦、波浪和漩涡3个部分进行分析，其中摩擦阻尼部分在总阻尼里占的比重最少，而漩涡阻尼部分占的比重最多，两者的相关力矩曲线具有一定的相位差。另外，波浪阻尼力矩的方向与摩擦和漩涡这2种阻尼力矩的方向是相反的。

本文还采用RANS方法对DTMB 5512三维船模进行了自由衰减的横摇运动模拟，验证了网格以及计算方法的可行性。在此基础上，对该三维裸船模与带舭龙骨船体进行了强迫横摇运动的模拟计算，将摩擦阻尼部分与舭龙骨阻尼部分从总阻尼中分离出来进行了分析，与二维结果类似，摩擦阻尼部分占比非常小，波浪与漩涡阻尼的影响较大。而舭龙骨则在一定的横摇角度范围内减小了总阻尼力矩，即阻尼力矩的舭龙骨影响部分曲线与总阻尼力矩曲线具有一定的相位差。

在以后的研究中，还将分析舭龙骨在船模自由横摇运动中的影响，并将欧拉方法应用于三维船模强迫横摇中，以分离得到波浪阻尼系数。同时，船模横摇相关参数的变化对于船模横摇阻尼系数各个部分的影响也有待进一步的研究。

［1］ CHAKRABARTI S.Empirical calculation of roll damp⁃ing for ships and barges［J］.Ocean Engineering，2001，28（7）：915-932.

［2］ IKEDA Y，HIMENO Y，TANAKA N.On eddy making component of roll damping force on naked hull［J］. Journal of the Society of Naval Architects of Japan，1977，142：54-64.

［3］ IKEDA Y，HIMENO Y，TANAKA N.Prediction meth⁃od for ship roll damping：00405［R］.Osaka：Depart⁃ment of Naval Architecture，University of Osaka，1978.

［4］ ALESSANDRINI B，DELHOMMEAU G.A fully cou⁃pled Navier-Stokes solver for calculation of turbulent incompressible free surface flow past a ship hull［J］.In⁃ternational Journal for Numerical Methods in Fluids，1999，29（2）：125-142.

［5］ LUCA B，STEFANO B.Influence of viscosity on radia⁃tion forces：a comparison between monohull，catamaran and SWATH［C］//Proceedings of the Twenty-third In⁃ternational Offshore and Polar Engineering.Anchorage， Alaska，USA：International Society of Offshore and Po⁃lar Engineers，2013.

［6］ MILLER R，GORSKI J，FRY D.Viscous roll predic⁃tions of a circular cylinder with bilge keels［C］//Pro⁃ceedings of the 24th Symposium on Naval Hydrodynam⁃ics.Fukuoka，Japan：NRNAS，2002.

［7］ WILSON R V，PATERSON E G，STERN F.Unsteady RANS CFD method for naval combatants in waves［C］// Proceedings of the 22nd ONR Symposium on Naval Hy⁃drodynamics.Washington，DC：［s.n.］，2006.

［8］ YANG C L，ZHU R C，MIAO G P，et al.Numerical sim⁃ulation of rolling for 3-D ship with forward speed and nonlinear damping analysis［J］.Journal of Hydrody⁃namics，Series B，2013，25（1）：148-155.

［9］ SHEN Z R，CAO H J，YE H X，et al.Manual of CFD solver naoe-FOAM-SJTU：2012SR118110［P］.Shang⁃hai：Shanghai Jiaotong University，2012.

［10］ JASAK H，JEMCOV A，ZELJKO T，et al.Open⁃FOAM：a C++library for complex physics simulations［C］//International Workshop on Coupled Methods in Numerical Dynamics.Dubrovnik.Croatia：IUC，2007：1-20.

［11］ ISSA R I.Solution of the implicitly discretised fluid flow equations by operator-splitting［J］.Journal of Computational Physics，1986，62（1）：40-65.

［12］ 吴建威，尹崇宏，万德成.基于三种方法的螺旋桨敞水性能数值预报［J］.水动力学研究与进展，2016，31（2）：177-187. WU J W，YIN C H，WAN D C.Numerical prediction of the propeller open-water performance based on three numerical methods［J］.Chinese Journal of Hy⁃drodynamics，2016，31（2）：177-187（in Chinese）.

［13］ SHEN Z R，WAN D C.An irregular wave generating approach based on naoe-FOAM-SJTU solver［J］.Chi⁃na Ocean Engineering，2016，30（2）：177-192.

［14］ 李鹏飞，万德成，刘建成.基于致动线模型的风力机尾流场数值模拟［J］.水动力学研究与进展，2016，31（2）：127-134.

LI P F，WAN D C，LIU J C.Numerical simulations of wake flows of wind turbine based on actuator line model［J］.Chinese Journal of Hydrodynamics，2016，31（2）：127-134（in Chinese）.

［15］ 尹崇宏，吴建威，万德成.基于IDDES方法的模型尺度和实尺度VLCC阻力预报与流场分析［J］.水动力学研究与进展，2016，31（3）：259-268.

YIN C H，WU J W，WAN D C.Model-and full-scale VLCC resistance prediction and flow field analysis based on IDDES method［J］.Chinese Journal of Hy⁃drodynamics，2016，31（3）：259-268（in Chinese）.

［16］ CAO H J，WAN D C.Development of multidirectional nonlinear numerical wave tank by naoe-FOAM-SJTU solver［J］.International Journal of Ocean System En⁃gineering，2014，4（1）：49-56.

［17］ ZHA R S，YE H X，SHEN Z R，et al.Numerical com⁃putations of resistance of high speed catamaran in calm water［J］.Journal of Hydrodynamics，2014，26（6）：930-938.

［18］ SHEN Z R，WAN D C，CARRICA P M.Dynamic over⁃set grids in OpenFOAM with application to KCS self-propulsion and maneuvering［J］.Ocean Engi⁃neering，2015，108：287-306.

［19］ 杨波，石爱国，吴明.基于计算流体力学理论的船舶横摇阻尼系数计算［J］.中国航海，2012，35（3）：76-80.

YANG B，SHI A G，WU M.Calculation of ship roll damping coefficient based on CFD［J］.Navigation of China，2012，35（3）：76-80（in Chinese）.

［20］ 李积德.船舶耐波性［M］.哈尔滨：哈尔滨工程大学出版社，2007.

LI J D.Ship seakeeping［M］.Harbin：Harbin Engi⁃neering University Press，2007（in Chinese）.

［21］ OLIVIERI A，PISTANI F，AVANZINI A，et al.Tow⁃ing tank experiments of resistance，sinkageand trim，boundary layer，wake，and free surface flow around a naval combatant INSEAN 2340 model：IIHR-TR-421［R］.Iowa：The University of Iowa，2001.

［22］ 黄昊，郭海强，朱仁传，等.粘性流中船舶横摇阻尼计算［J］.船舶力学，2008，12（4）：568-573.

HUANG H，GUO H Q，ZHU R C，et al.Computations of ship roll damping in viscous flow［J］.Journal of Ship Mechanics，2008，12（4）：568-573（in Chinese）.

［23］ KORPUS R A，FALZARANO J M.Presiction of vis⁃cous ship roll damping by unsteady Navier-Stokes techniques［J］.Journal of Offshore Mechanics and Arctic Engineering，1997，119：108-113.

Numerical analysis of viscous effect on ship rolling motions based on CFD

LUO Tian1，2，3，WAN Decheng1，2，3
1 School of Naval Architecture，Ocean and Civil Engineering，Shanghai Jiao Tong University，Shanghai 200240，China
2 State Key Laboratory of Ocean Engineering，Shanghai Jiao Tong University，Shanghai 200240，China
3 Collaborative Innovation Center for Advanced Ship and Deep-Sea Exploration，Shanghai 200240，China

During the ship design procedure,the analysis of ship rolling motions is of great significance because the rolling motions have extraordinary effects on the sea-keeping,maneuverability and stability of a ship.It is difficult to simulate rolling motions due to the effect of viscosity,which causes many nonlinear components in computation.As such,the potential theory used for other ship motions cannot be used for rolling motions.This paper simulates the rolling motions of the DTMB 5512 ship model and the ship transverse section of the S60 ship model with a naoe-FOAM-SJTU solver using the Reynolds Averaged Navier Stokes（RANS）method based on the OpenFOAM.The results of rolling motions are compared with the experimental data,which confirms the reliability of the meshes and results.For the ship transverse section of the S60 ship model,the damping coefficient is divided into three parts with the Euler and RANS methods:friction,vorticity and wave parts.For the DTMB 5512 ship model,the damping coefficient is also respectively analyzed,including the friction,vorticity,wave and bilge keel parts.The results in this paper show that the vorticity part accounts for the greatest proportion,while the friction part accounts for the least, and the bilge keels reduces the damping moment to a certain extent which shows the effect of rolling parameters on rolling motions and moments.

ship；rolling motions；damping coefficient；damping moment；RANS method；Euler method

U661.3

A

10.3969/j.issn.1673-3185.2017.02.001

http://kns.cnki.net/kcms/detail/42.1755.TJ.20170313.1618.036.html

罗天，万德成.基于CFD的船舶横摇数值模拟与粘性效应分析［J］.中国舰船研究，2017，12（2）：1-11，48. LUO T，WAN D C.Numerical analysis of viscous effect on ship rolling motions based on CFD［J］.Chinese Journal of Ship Research，2017，12（2）：1-11，48.

2016-09-01 < class="emphasis_bold"> 网络出版时间：

时间：2017-3-13 16:18

国家自然科学基金资助项目（51379125，51490675，11432009，51579145）；长江学者奖励计划（T2014099）；

上海高校东方学者特聘教授岗位跟踪计划（2013022）；工信部数值水池创新专项VIV/VIM项目（2016-23/09）

罗天，女，1991年生，硕士生。研究方向：粘性流场中船舶的横摇运动分析。

E-mail：sky360652272@163.com

万德成（通信作者），男，1967年生，博士，教授。研究方向：船舶与海洋工程水动力学及计算流体力学的基础理论及其应用。E-mail：dcwan@sjtu.edu.cn

期刊网址：www.ship-research.com