0 引　言

船舶操纵性直接关系到船舶的安全航行、高效作业以及环境保护等多方面因素，船舶操纵性预报是船舶设计过程中重要的一个环节，随着计算机水平的提高，数值模拟技术的应用愈发成熟，船舶操纵性数值模拟研究也成为近年来船舶领域极为重要的研究方向。数值模拟技术能够在设计阶段就预测和改进船舶的操纵性能，这对于降低实船试验成本、缩短设计周期具有显著效果。高精度的数值模拟，能够在复杂海况下对船舶操纵性进行深入分析，为船舶避免碰撞、减少事故提供科学依据，具有实际的意义。

近年来，大量的国内外学者们开展了真实的自航船模操纵性数值模拟研究，也取得了一些相应的研究成果。学者们在处理船-桨-舵耦合的方面对螺旋桨的模拟形式不同，有的学者采用以恒定推力代替螺旋桨以简化计算来使得船舶自航，有的学者采用螺旋桨体积力的方法来替代真实的螺旋桨建模，在条件更充足的情况下，也有部分学者对螺旋桨真实建模，实现船-桨-舵的全耦合自航操纵预报数值模拟。

操戈等^[1]利用自研URANS求解器HUST-Ship对船模ONRT在波浪中的回转操纵运动进行直接CFD仿真，通过与Tokyo2015和Simman2020提供的数据进行对比，可以得到自研CFD软件对操纵性的求解很适用，而且迭代型体积力法较传统体积力法更好地反映螺旋桨处的流场信息。王子平^[2]利用重叠网格方法在STAR-CCM+平台对实尺度的DTMB-5514船舶操纵性预报，对真实的桨和舵也进行建模，预报出其在回转试验和Z形操纵试验过程中的六自由度运动时历，经过对比该方法与试验数据吻合良好。朱锋等^[3]在自研粘流计算软件Oship上采用体积力法对高速巡逻艇的回转性能进行了预报，同时考虑了扰流板对其回转性能的影响。王建华^[4]基于上海交通大学开发的naoe-FOAM-SJTU水动力求解器，又加入了船舶自航操纵控制模块，实现了船、桨、舵相互作用下的船舶操纵运动数值模拟。郑秋萌^[5]基于CFD数值模拟的方法，利用螺旋桨体积力模型模拟了自航船模KVLCC2的静水回转运动，对其六自由度运动和纵向以及横向水动力值进行预报，并为其后续神经网络预报打下基础。于东等^[6]利用恒定的推力和力矩来模拟一艘三体船的推进器，在STAR-CCM+平台模拟了该三体船的回转试验，分析了其回转性能。

Choi等^[7−8]在商业软件Fluent平台结合自开发的螺旋桨体积力模型，模拟了多种带螺旋桨的船舶的自航推进性能。Simonsen等^[9]通过结合势流理论与体积力模型，研究了螺旋桨对船舶作用的数值模拟方法，该研究采用迭代技术，实现了RANS方法计算得到的流场与螺旋桨产生的体积力之间的耦合，此外通过该方法模拟了带附体的油船Esso Osaka的操纵运动，详细分析了船—桨—舵间的相互作用效应。Sukas等^[10]采用直接CFD模拟方法，基于STAR-CCM+平台的螺旋桨体积力模型，直接预报了DTC船模的回转和Z形操纵运动的六自由度时历，并且和系统建模的方法对比，在网格数量小的条件下，直接CFD的方法时间和精度上都取得了不错的效果，网格量大的时候，CFD方法时间成本增大，但是能精确给出流场数据。Kim等^[11]也基于体积力法和重叠网格方法模拟了KVLCC2船在波浪下的回转实验，模拟发现，短波对回转能力的影响最为显著，波浪会影响转弯性能，但是稳定回转的回转特性变化不大。

综上所述，学者们在船舶自航操纵数值模拟方面已经开展了大量研究，大部分研究都是以体积力法来代替螺旋桨的建模，但是大多数研究都是对国际标模的研究，没有针对实际的工程应用问题，而且只有少部分的学者考虑了波浪条件下的操纵性数值模拟分析与静水条件下的对比，本文为了解决某双船首构型的风电运维母船的操作性能预报并为后续优化设计打下基础，基于体积力法和重叠网格方法展开回转试验研究，相关研究成果可为新型风电运维母船设计提供工程意义。

1 操纵性计算CFD数值方法 1.1 首尾对称风电运维母船模型

本文研究对象为首尾结构对称的风电运维母船，即双船首的构型，其模型如图1所示，该船模与实船比例为1∶25，其主尺度表如表1所示，缩尺后的模型的Lpp=3.152 m，本文全部计算都使用此缩尺模型，为了模拟真实自航状态下的船舶回转试验，需要对推进器进行模拟，本船采用传统的双桨双舵推进方式，对舵进行真实建模，下文将介绍控制其旋转的方法以模拟真实的回转试验，而螺旋桨采用体积力法模型替代，本文计算全部在STAR-CCM+平台完成，选用STAR-CCM+平台的虚拟

1.2 数值模拟基本方法 1.2.1 控制方程

流体的基本控制方程由连续性方程、动量方程以及能量方程构成，它们分别对应于质量守恒定律、牛顿第二运动定律和能量守恒定律。由于本文的数值模拟不涉及能量的转换，这里不考虑能量方程。

1）连续性方程

连续性方程源自于质量守恒原理。在流场中任取任意形状的任意空间作为控制体，该控制体的体表面称为控制面，流场中的流体通过该表面流入控制体，也通过该表面流出。

根据质量守恒定律有，控制体内部流体质量的变化量等于单位时间内流入和流出该控制体的流体质量的差值。由此可以得到流体的连续性方程表达式为：

$ \frac{\partial }{{\partial t}}\iiint\limits_V {\rho {\mathrm{d}}x{\mathrm{d}}y{\mathrm{d}}z} + \mathop{\oiint}\limits_A {\rho v \cdot n{\mathrm{d}}A} = 0。$

(1)

式中：$V$为控制体；$A$为控制面。

在笛卡尔坐标系中，式(1)可表示为：

$ \frac{\partial p}{\partial t}+u\frac{\partial(\rho u)}{\partial x}+v\frac{\partial(\rho v)}{\partial y}+w\frac{\partial(\rho w)}{\partial z}=0。$

(2)

式中：$u$、$v$、$w$分别为速度矢量$\overrightarrow V $在笛卡尔三维坐标系下的3个不同坐标轴方向的速度分量。

当流体为恒密度均匀不可压流体时，该方程可以简化为：

$ \frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}+\frac{\partial w}{\partial z}=0。$

(3)

在柱坐标系中，其形式为：

$ \frac{{{v_r}}}{r} + \frac{{\partial {v_r}}}{{\partial r}} + \frac{{\partial {v_\theta }}}{{r\partial \theta }}{\text{ + }}\frac{{\partial {v_z}}}{{\partial z}}{\text{ = 0}} 。$

(4)

2）动量方程

流体力学领域基于牛顿第二定律推导出来的方程为动量方程。纳维-斯托克斯方程（Navier-Stokes equations），通常简称为N-S方程，其定义为：流体微元的动量时间变化率等于作用在该流体微元上的力，对应的表达式为：

$ \left\{ \begin{gathered} \rho \frac{{{\mathrm{d}}u}}{{{\mathrm{d}}t}} = \rho {F_{bx}} + \frac{{\partial {p_{xx}}}}{{\partial x}} + \frac{{\partial {p_{yx}}}}{{\partial y}} + \frac{{\partial {p_{zx}}}}{{\partial z}}, \\ \rho \frac{{{\mathrm{d}}v}}{{{\mathrm{d}}t}} = \rho {F_{by}} + \frac{{\partial {p_{xy}}}}{{\partial x}} + \frac{{\partial {p_{yy}}}}{{\partial y}} + \frac{{\partial {p_{zy}}}}{{\partial z}}, \\ \rho \frac{{{\mathrm{d}}z}}{{{\mathrm{d}}t}} = \rho {F_{bz}} + \frac{{\partial {p_{xz}}}}{{\partial x}} + \frac{{\partial {p_{yz}}}}{{\partial y}} + \frac{{\partial {p_{zz}}}}{{\partial z}}。\\ \end{gathered} \right. $

(5)

式中：${F_{bx}}$、${F_{by}}$、${F_{bz}}$是质量力在笛卡尔立体坐标系内的3个分量。

一般在船舶流体力学的研究中，所研究的流体的质量密度和粘性系数都为常数，此时的N-S方程为：

$ \rho\frac{\mathrm{d}v}{\mathrm{d}t}=\rho F-\rm{grad}\rho+\mu\nabla^2v。$

(6)

对于理想流体来说，忽略流体的粘性，则可得其动量方程为：

$ \left\{ \begin{gathered} \frac{{{\mathrm{d}}u}}{{{\mathrm{d}}t}} = \frac{{\partial u}}{{\partial t}} + u\frac{{\partial u}}{{\partial y}} + w\frac{{\partial u}}{{\partial z}} = {f_x} - \frac{{\partial p}}{{\rho \partial x}}, \\ \frac{{{\mathrm{d}}y}}{{{\mathrm{d}}x}} = \frac{{\partial v}}{{\partial t}} + u\frac{{\partial v}}{{\partial y}} + w\frac{{\partial v}}{{\partial z}} = {f_y} - \frac{{\partial p}}{{\rho \partial y}} 。\\ \end{gathered} \right. $

(7)

1.2.2 湍流模型

当今主要的湍流数值模拟方法可以分为直接数值模拟（DNS）和非直接数值模拟2种，其中非直接数值模拟方法主要由大涡模拟法（LES）和雷诺时均法（RANS）方法构成。DNS方法能够无需依赖于特定的湍流模型，通过直接求解时间相关的N-S方程来精确模拟湍流的所有尺度，从而提供详尽的湍流结构信息；相较之下，LES方法对湍流中的大尺度涡流进行直接模拟，而小尺度涡流则通过引入封闭模型来处理，通过应用滤波函数实现对不同尺度涡流的分离，旨在简化求解过程而保留关键的湍流特性。然而目前在船舶CFD研究领域，应用直接模拟法和大涡模拟法的较少。

雷诺时均法其特点在于求解的效率很高，精度经过许多学者的验证也满足许多工程领域研究的需求，是目前船舶水动力研究领域最常用的湍流数值模拟方法。其主要思想是求解时均化的雷诺（Reynolds）方程，通过模型在时均化的方程表示瞬态的脉冲量。Reynolds方程表达式为：

$ \left( {\rho {u_i}} \right) + \frac{\partial }{{\partial {x_j}}}\left( {\rho {u_i}{u_j}} \right) = - \frac{{\partial p}}{{\partial {x_j}}} + \frac{\partial }{{\partial {x_j}}}\left( {\mu \frac{{\partial {u_i}}}{{\partial {x_j}}} - \rho \overline {u_i^{'}u_j^{'}} } \right) + {S _i} 。$

(8)

在雷诺时均法中，按对雷诺应力的处理方法的不同，主要可以分为雷诺应力模型和涡粘模型两大类。涡粘模型引入了Boussinesq提出的涡粘假设：湍流脉动所造成的附加应力也与层流运动应力那样可以同时均的应变率关联起来，该涡粘模型表达式为：

$ - \rho \overline {u_i^{'}u_j^{'}} = {\mu _t}\left( {\frac{{\partial {u_i}}}{{\partial {x_j}}} + \frac{{\partial {u_j}}}{{\partial {x_i}}}} \right) - \frac{2}{3}\left( {\rho k + {\mu _t}\frac{{\partial {u_i}}}{{\partial {x_i}}}} \right){\delta _{ij}}。$

(9)

式中：$ {\mu }_{t} $为湍流粘度；$ {\mu }_{i} $为时均速度，当$ i=j $时，$ {\delta }_{ij}=1 $；当$ i\ne j $时，$ {\delta }_{ij}=0 $；$ k $为湍动能。

加入Boussinesq假设后，确定湍流粘度$ {\mu }_{t} $成为模拟计算湍流的关键问题，在应用中最常见的是两方程式涡粘模型包括标准$ k-\omega $模型、标准k − ε模型、RNG k − ε模型以及Realizable k − ε模型等，本文计算主要采用两方程SST $ k-\omega $模型。

1.3 自航推进和回转控制模块 1.3.1 螺旋桨体积力模型

本文采用STAR-CCM+的虚拟盘模型的体积力方法模拟螺旋桨对船体的推力，利用体积力法，优点是可以准备模拟螺旋桨对船体的作用力并能减小网格数量，减低计算成本，缺点是该方法无法获得螺旋桨周围的详细流场，但是本文目的在于研究船舶的操纵性，对螺旋桨周围流场的详细信息不予考虑，为了节约计算资源使用虚拟盘体对船舶产生自航推进的作用。

体积力螺旋桨法对螺旋桨的推力和扭矩进行建模，从而在不实际求解螺旋桨几何的情况下模拟推进力。该方法将体积力$ {f}_{b} $均匀分布在虚拟盘体上，体积力沿径向分布：

$ \left\{ {\begin{array}{*{20}{l}} {{f_{bx}} = {A_x}{r^*}\sqrt {1 - {r^*}} } ，\\ {{f_{b\theta }} = {A_\theta }\dfrac{{{r^*}\sqrt {1 - {r^*}} }}{{{r^*}\left( {1 - r_h^{'}} \right) + r_h^{'}}}} 。\end{array}} \right. $

(10)

式中：$ {f}_{bx} $为轴向体积分力，N；$ {f}_{b\theta } $为切向体积分力，N；$ {r}^{\mathrm{*}} $和$ {r}'_{h} $为无量纲化的盘径，可以表示为桨盘轮毂半径和叶端半径的函数，如下式所示：

$ \left\{ {\begin{array}{*{20}{l}} {{r^*} = \dfrac{{r' - r_h^{'}}}{{1 - r_h^{'}}}}, \\ {r_h^{'} = \dfrac{{{R_H}}}{{{R_P}}}} ,\\ {r' = \dfrac{r}{{{R_p}}}} 。\end{array}} \right. $

(11)

式中：r为径向坐标，m；$ {R}_{H} $为轮毂半径，m；$ {R}_{p} $为螺旋桨叶尖端半径，m。

常数$ {A}_{x} $和$ {A}_{\theta } $用下式计算：

$ \left\{ \begin{gathered} A_{x}=\dfrac{105}{8} \dfrac{T}{{\text{π}} \Delta\left(3 R_{H}+4 R_{P}\right)\left(R_{P}-R_{H}\right)}，\\ A_{\theta}=\dfrac{105}{8} \dfrac{Q}{{\text{π}} \Delta\left(3 R_{H}+4 R_{P}\right)\left(R_{P}-R_{H}\right)}。\end{gathered}\right. $

(12)

式中：$ \Delta $为虚拟盘体厚度，m；T为推力，N；Q为扭矩，N·m。

因此，虚拟盘体积力模型需要确定盘的尺寸和方向以及螺旋桨敞水特性曲线，同时需要将虚拟盘体局部坐标系关联在船体随体坐标系下以实现对船体的推力作用。

1.3.2 船舵控制方程

STAR-CCM+平台可以利用多级物体运动模块，将船舶的体运动叠加一个坐标系嵌套的舵旋转运动，以控制舵的运动改变船舶的运动状态，同时为了确保物理量能够在不同计算域之间传递，需要将舵区域与重叠域以及背景域建立多个重叠网格交界面。

船舶回转运动舵的旋转将被设置为有关于时间和旋转速率的函数，回转运动需要实现舵角从0°转到最大舵角的控制，其控制方程为：

$ {\delta }^{{'}}(t)=\left\{ {\begin{array}{*{20}{l}}0\text{，}{t}_{0} > t > 0，\\ \pm \Delta {\delta }_{\mathrm{max}}\text{，}{t}_{1} > t > {t}_{0}，\\ 0\text{，}{\mathit{t}}_{2} > t > {t}_{1}。\end{array}} \right. $

(13)

式中：${t_0}$为船舶直航的时间；${t_1}$为船舶操舵完成时间；${t_2}$为船舶回转试验结束时刻；$\Delta {\delta _{\max }}$为最大转舵速率。

1.4 数值方法验证

本文选用国际标模DTMB 5415船模的纯横荡运动数值模拟来验证本文的数值模拟方法。船舶纯横荡运动即在纵向运动的同时叠加一个横向的正弦振荡运动，属于动态约束船模试验。使用的船模缩尺比为1∶24.830，船模主要参数为$ Lpp=5.719 $ m，$ B= 0.768 $ m，$ T=0.248 $ m，$ \mathrm{\Delta }=0.554 $ m³。其中，$ Lpp $、$ B $和$ T $分别为船模垂线间长、水线宽和吃水，$ \mathrm{\Delta } $为排水量。在STAR-CCM+平台划分计算域并采用本文的数值模型进行计算，计算工况为航速$ V=2.097 $ m/s，横荡频率$ f=0.098 $ Hz，横荡振幅$ a=0.594 $ m。

计算对船体的横向力进行监测，将数据导出后首先利用傅里叶级数进行拟合，可以得到相应的水动力系数，并求得水动力导数，将这些水动力导数值与simman2008国际水会所提供的试验值进行对比，如表2所示。从结果对比来看，水动力导数的误差较小，然后将前10 s的数据与试验数据进行对比，见图2。可以看出总体来说数据趋势一致，尽管前几秒由于流场不稳定导致极值有所偏差，计算结果可以接受，说明本文所采用的数值模拟方法在动态船模试验中相当可靠。

2 自航运动数值模拟 2.1 工况设置与求解设置

首先分析该船舶在自航状态下的六自由度运动情况，还要在此基础上求出船舶在最大服务航速下的自航点，以便后续开展船舶回转操纵试验数值模拟。根据螺旋桨转速的不同，设置了3组不同的自航工况，如表3所示，螺旋桨转速分别为5、7.5、10 r/s。该船舶的最大服务航速为14.5 kn，对应的傅汝德数为0.267。

自航点的求解时需要模拟船舶在静水域中的完全自航状态，需要使用自航推进器来使得船舶向前推进，同时依靠船舵进行航向保持，所以计算域设置时需要考虑到舵的区域和推进器区域，而此处推进器使用螺旋桨体积力法即STAR-CCM+平台的虚拟盘体代替，所以不需要设置专门的计算域，因此全部的计算域设置为：试验水池区域即背景网格区域、船体周围区域即重叠网格区域以及左舵区域和右舵区域。背景区域与船体周围区域生成重叠网格，左舵区域分别与背景区域和船体周围区域生成重叠网格，右舵区域分别与背景区域和船体周围区域生成重叠网格。其中背景域大小可设置为如图3所示，在船长方向，设置船尾后$ 2.5{L}_{pp} $水域，船首前$ 1.5{L}_{pp} $水域；在船侧区域分别设置$ 2{L}_{pp} $大小，在船深方向，设置吃水以下$ 1{L}_{pp} $，吃水以上$ 0.5{L}_{pp} $。船体周围区域设置如下：船长方向前后各$ 0.2{L}_{pp} $，船宽方向距离边缘0.3B，高度方向也为0.3B。流域的入口和底部均设置为速度入口条件，顶部和出口设置为压力出口边界条件，侧面均设为对称平面的边界条件。左右舵区域只需要保证有3～5层过渡网格大小区域即可，不做过多设计。

本文设置了5个重叠网格交界面，背景域网格和船周重叠网格区域设置基础尺寸为$ 0.02{L}_{pp} $，需要对自由液面处和开尔文尾流区以及重叠网格过渡区域进行加密设置，自由液面区域即在吃水深度处进行加密即可，因为没有波浪不需要过于密，节约计算资源，同时开尔文尾流区也尽量小，为了保证重叠区域之间物理量的交换不出错，在交界区域设置5层过渡网格，网格划分结果总数为365万，其中船体周围重叠域250万，水池背景域和其他区域115万，网格整体如图4所示。

2.2 计算结果分析

在计算过程中分别监测虚拟盘体产生的推力和船体所产生的阻力，它们二者的差值就是船舶在纵向受到的净力值，船模的自航点就是在设定航速下达到动平衡的状态，根据牛顿第二定律，此时的净力值应该为0。表3中的3组工况船模初始航速为V₀=1.49 m/s（对应实船14.5 kn），监测到的净力值分别为−1.915 、−0.632 、1.384 N。根据该数据绘制净力与转速的曲线，如图5(a)所示，在图中拟合曲线与X轴的交点，此处对应的螺旋桨转速即为自航点，即约为n=8.35 r/s。图5(b)为船体受到净力的结果，船体所受的净力大小基本在平衡点附近波动，说明此时船舶螺旋桨产生的推力和此航速下的船舶所受的阻力基本相等，船舶达到动平衡状态。

3 回转运动数值模拟 3.1 工况设置与求解设置

本文使用数值模拟的方法来模拟该船舶的回转试验，回转试验也同样采用L3的缩小比例的模型来模拟，为了全面分析这种结构的船舶的回转性能，设置了静水和波浪条件下的计算，来分析该船舶在静水和波浪条件下回转的差异性。如表4所示共设置了6个工况，静水左右35°舵角、迎浪35°舵角以及随浪35°舵角，初始航速为V₀=1.49 m/s，取左舵为舵角负值，波长$ \lambda ={L}_{WL} $，波高H/$ \lambda $=0.02，记该波浪工况为wave1。

回转运动的数值模拟，在自航运动的基础上，需要控制舵的旋转，在STAR-CCM+平台根据上文所提到的舵控制的方程，编写UDF来控制舵的旋转，利用多级物体运动模块从而实现船舶的回转。静水条件下的计算域设置可以和静水条件下的自航完全相同，对背景域的控制，需要控制其与船舶保持相对地运动，以此来节约计算资源，即控制背景域跟随船体回转。而波浪条件下的回转运动，计算域大小保持与静水条件下相同，但是网格需要对波高范围内加密，对兴波面进行更细化的加密，网格加密区如图6所示，最终全域的网格总数在590万左右。

3.2 计算结果分析

首先对该首尾对称型船舶的静水回转轨迹进行分析，计算过程中监测到的回转轨迹如图7所示，为了对比分析一般船舶与这艘船舶的回转性能的差异，对回转轨迹中的每个时刻的轨迹点都无因次化处理。这艘船的回转轨迹可以看出，在初始回转时刻其转首性能较差，其战术回转直径较一般船舶较大，而稳定回转半径在工程船舶中算极好的情况，是由于其首尾对称的结构导致在回转过程中受到尾流场的影响较小，通过分析得到，左右舵回转的定常回转半径不完全相同，本船左舵的回转要略差于右舵，这与许多实船试验的结果也类似。静水回转参数和波浪下的回转参数如表5所示，由于这艘船舶为肥大型船舶，其转首性能较差，其纵距和正横距较大，即初始转首性能较差。

在波浪条件下，由于受到波浪的影响，其回转性能也受到了影响，从图8中可知，在随浪的条件下，回转半径会比静水略小一些，而迎浪的条件和静水中回转半径相差不大几乎一致，但随浪导致纵距增大了约6%的船长，迎浪条件下减小了约5%的船长，这是由于波浪分别对回转运动的初速度衰减有一定的减小或增加。综上可知，该船舶在静水和波浪的条件下，都表现出了良好的回转性能，但是其缺点在于转首性能较差。

4 结　语

本文在STAR-CCM+平台基于体积力方法和重叠网格方法对一艘首尾结构对称的风电运维母船的缩尺模型进行回转操纵试验的研究，得到结论如下：

1）本文采用的数值模拟方法能够很好地模拟船舶自航操纵运动，通过模拟得出本文所研究的双船首构型的风电运维母船其回转性能较一般工程船舶和运输船舶优越。

2）在随浪和迎浪的条件下，该船舶的稳定回转性能受到的干扰较小，可见其在波浪中也具有良好的回转操作性能。

3）该船舶同样具有着肥大型船舶的缺点，其初始转首性能较差，战术回转直径较大，若要优化其操纵性能可着重优化其转首性能，如增加首推进器等，则可以得到操纵性能全面优秀的船型。