0 引　言

根据海运经济与物流研究所（ISL）数据，2024年全球集装箱吞吐量指数攀升至140.03，总运输量达5312万TEU，双双创下历史新高^[1]。这表明国际集装箱航运市场仍旧保持着上升状态，大型集装箱船运输依旧是航运市场的中坚力量。随着深圳西部港区集装箱运输的快速发展，船舶大型化趋势对港口基础设施提出了更高要求。为满足大型集装箱船全天候安全通航需求，不仅需要深水泊位，还必须保证足够的富裕水深（Under Keel Clearance, UKC），同时粤港澳大湾区建设和深圳先行示范区政策实施后，深圳港超大型船舶吞吐量显著增长。又受香港青马大桥净空限制^[2]，基于此铜鼓航道凭借其优越的通航条件，正日益成为大型集装箱船进出深圳西部港区的首选航线之一。因此，为了满足大型集装箱船的通航需求，铜鼓航道向两侧扩建，这也导致航道两侧岸壁倾斜，同时航道内的风流作用也使得大型集装箱船的通航难度大大增加。

计算流体力学（Computational Fuild Dynamics，CFD）凭借计算机技术和数值方法的发展，已成为研究船舶在受限水域运动规律的重要工具。张金鹏^[3]通过对比多种风载荷系数计算方法对集装箱船的应用效果，验证了数值模拟在计算集装箱船风载荷方面的有效性。谢海波等^[4]基于三维势流理论，采用SESAM软件计算了集装箱船、油船和散货船在不同波浪条件和水深吃水比下的锚泊下沉量，但只依赖势流理论，未考虑粘性效应和湍流影响。陈明达等^[5]基于雷诺平均纳维-斯托克斯模型（RANS）方法和动网格技术，验证了CFD方法在限制水域水动力计算中的可行性，揭示了航道形状不对称性对船舶受力的显著影响，但未考虑船舶转向或斜航等更复杂的操纵状态。GUO等^[6]虽然通过CFD与物理模型试验探讨了阻塞效应，但仅模拟船舶与岸壁距离较大的开放水域，忽略了浅水效应与岸壁效应的耦合作用。CHILLCCE 等^[7]采用三种流场（无黏流、层流、湍流）对比分析不同水深下黏性力对船舶纵倾和下沉的影响，明确黏性效应显著影响纵倾，湍流模型（URANS）在浅水中计算更准确。Lee等^[8]使用CFD方法分析了2种不同大型船舶近岸航行的水动力，并研究了各种因素对水动力相互作用的影响。Hadi等^[9]发现航道约束程度增加会使船体周围水流压力分布发生变化，导致船舶下沉量增大。夏明宇等^[10]基CFD方法分析了海岬型散货船在不同流速下流致船体下沉量。毛超进等^[11]使用CFD方法探究不同船速、不同水深吃水比以及不同船岸距离下对船舶流压场变化情况。综上，CFD技术能够精确模拟复杂流场特性（如黏性效应、湍流分离、自由液面变形等），弥补势流理论忽略粘性和涡旋的缺陷，但主要针对船舶直航不受斜风斜流的状态，对大型集装箱船在限制水域（铜鼓航道特点的水域）研究较少。

本文基于黏流理论，利用 SolidWorks 软件对铜鼓航道进行 3D 建模；在杜伊斯堡测试案例集装箱船（Duisburg Test Case, DTC）的基础上，通过 SolidWorks 对甲板集装箱、上层建筑进行3D建模；应用STAR-CCM+（URANS）求解器，对DTC在的铜鼓航道航行受斜风、斜流致船舶偏转、下沉等进行数值模拟；从而得到DTC在铜鼓航道典型水文气象条件下航行时船舶运动姿态的变化规律，旨在为大型集装箱船在铜鼓航道内安全航行提供参考。

1 数值模拟方法及研究对象 1.1 控制方程

雷诺数（Reynolds Number）用来表示流体流动情况的无量纲数，流体根据雷诺数的大小分为层流和湍流两种流动状态，其中湍流对船舶的运动会产生明显的影响作用。为准确模拟黏性流体及湍流效应，本研究采用雷诺平均方法（Reynolds-averaged Approach），通过时间平均与涡度分解对流动过程进行建模。控制方程选用三维不可压缩流体的雷诺平均纳维-斯托克斯方程（Reynolds-averaged Navier-Stokes，RANS），该方程组完整刻画了流体的连续性、动量传递及能量输运特性。方程组如下所示^[12]：

$ \frac{\partial \left(\rho \overline{{U}_{i}}\right)}{\partial {x}_{i}}=0，$

(1)

$ \frac{\partial \left(\rho \overline{{u}_{i}}\right)}{\partial t}+\frac{\partial }{\partial {x}_{j}}\left(\rho \overline{{u}_{i}}{\overline{u}}_{j}+\rho \overline{u_{i}^{'}}\overline{u_{j}^{'}}\right)=-\frac{\partial p}{\partial {x}_{i}}+\frac{\partial {\tau }_{ij}}{\partial {x}_{j}}，$

(2)

$ {\tau }_{ij}=\mu \left(\frac{\partial \overline{{u}_{i}}}{\partial {x}_{j}}+\frac{\partial {\overline{u}}_{j}}{\partial {x}_{i}}\right)。$

(3)

式中：$ {x}_{i\left(j\right)}(i,j=1{,}2,3) $为三维空间向量的各分量指标；$ {\overline{u}}_{i\left(j\right)} $为速度场的时均值；$ P $为时均压力场；$ \rho $为流体密度特性；$ \overline{u_{i}^{'}u_{j}^{'}} $为湍流脉动相关的雷诺应力张量；$ {\tau }_{ij} $为应力张量的各分量；μ为流体的动力粘性系数。

为了准确捕捉DTC在铜鼓航道进港航行时航道水域表面和空气界面的交界面的变化，引入VOF（Volume of Fluid）模型^[13]追踪和定位液-气2种流体之间的自由表面。VOF模型采用欧拉方法，在固定网格中求解动量方程组，并通过追踪各流体的体积分数来模拟不混溶流体间的界面。此方法引入一个标量函数——体积分数$ {\alpha }_{m} $，用于描述流体在控制体积（即计算网格单元）中的占比情况：定义$ {\alpha }_{1} $、$ {\alpha }_{2} $为两种流体1、2，若单元格完全为空，则$ {\alpha }_{1}=0 $；若单元格被流体完全占据，则$ {\alpha }_{2}=1 $；若单元格内存在流体界面，在原方程组中引入体积连续性方程使网格单元内存在两种流体的交界面即0 < $ {\alpha }_{2} $ < 1的区域。为保证液-气两种流体之间的自由表面即自由液面的准确性，使用CFD软件STAR-CCM+的高分辨率界面捕获（HRIC）。

$ \frac{\partial {\alpha }_{m}}{\partial t}+\frac{\partial \left({\alpha }_{m}{u}_{i}\right)}{\partial {x}_{i}}=0,\sum \limits_{m=1}^{2}{\alpha }_{m}=1。$

(4)

本文中动量方程中的未封闭项使用剪切应力输运SST k-ω湍流模型进行处理，该模型作为低雷诺数（$ Re $<200）的RANS封闭方案，能够精确捕捉近壁面区的流体动力特性^[14]。在数值求解过程中，构建分离流模型以非耦合形式求解流动方程，其中压力-速度耦合通过SIMPLE算法实现隐式迭代，并通过压力校正方程确保速度满足质量守恒的约束。

1.2 DTC计算模型

本文采用的船模为船模缩尺比为1∶89.11的杜伊斯堡测试案例集装箱船（Duisburg Test Case, DTC），由美国船舶技术、海洋工程和运输系统研究所（ISMT）在杜伊斯堡开发^[15]，该案例具备丰富的实验和模拟数据，用于进行比较和验证。本文采用满载的DTC船舶、选择拖曳方式进行数值模拟，因此忽略了螺旋桨与船舵的建模影响，几何形状如图1(a)，为了更好地去模拟风对集装箱船箱体与船体的载荷作用，通过SolidWorks对甲板集装箱、上层建筑进行3D建模，如图2(a)，建模符合船舶构造标准。所用船型比例与弗兰德斯水力学研究所采取的船型相同。并将模拟结果与弗兰德斯水力学研究所（FHR）进行物模实验的结果进行对比验证，以此来验证本文所采用数值方法的准确性和可行性。DTC船的详细信息见表1。

2 航道水域计算域 2.1 铜鼓航道水域环境计算域建模

铜鼓航道段长约 16 km，总体呈 NNE-SSW 走向，铜鼓航道设计底宽270 m，设计底标高−17.5 m，航道边坡为1∶7。基于铜鼓航道概况，为了方便计算，航道长度截取DTC实船的5倍船长，航道宽度及边坡高度按照DTC船模的缩尺比89.11换算后得到铜鼓航道模型尺寸，使用SolidWorks对其缩尺后的尺寸进行3 D建模。在航道内水域，湍流强度采用低值进行计算为1%，它用于表示速度脉动量的均方根与平均速度之间的比值。

图2为DTC在铜鼓航道航行计算模型示意图，其中入口边界、航道两侧、顶部为速度入口，出口边界为压力出口，船表面和航道底部及两侧斜壁为壁面。根据铜鼓航道目前的通航制度：3万吨级及以上船舶通过铜鼓航道实行单向通行^[16]，因此，计算域中，将船放置于航道水域中间。本次采用平波速度来模拟水流的速度和方向，入口边界距船中横剖面1.5$ L_{pp} $，出口边界距离船中横剖面2.5$ L_{pp} $。分别在入口、出口以及侧边进行长度为2.4 m的阻尼，确保在边界处的波浪更加接近实际情况，避免因边界反射造成的数值误差。

2.2 网格划分与时间步长

数值模拟的精度和收敛性高度依赖于网格质量。本研究采用STAR-CCM+中的切割体网格生成器，构建了非结构化六面体网格体系对计算域进行离散化处理。以h/T=1.2工况（见图3）为例，展示了船体纵剖面及自由液面区域的网格分布，其中对船首、船尾、底部边界层和自由液面进行网格加密。为提高计算效率与精度，增加多个控制体动态调节网格过渡区的密度梯度，使网格过渡更加合理，以此确保船舶湿表面和船首、尾波浪区域的流场特征（如边界层分离、自由面波动）能被精确捕捉。

数值模拟中，湍流的流动受壁面效应的影响较大，为更好的捕捉船体与航道底部及岸壁之间的水流的流动特征，建立壁面函数，引入一个无量纲参数$ {y}^{+} $，调整第1层边界层高度使近壁面$ {y}^{+} $≈1，设置1.5倍延伸比的6层边界层网格，$ {y}^{+} $的定义为^[17]：

$ {y}^+=y\mu \slash \left({U}_{\tau }\rho \right)。$

(5)

式中：$ y $为第一层边界层节点高度；$ \mu $为流体动力黏度；$ {U}_{\tau } $为壁面剪切应力估计速度；$ \rho $ 为流体密度。

UKC不足可能会导致船底与航道底部之间的边界层发生碰撞，从而产生更大的误差，因此对于浅水域问题，遵守CFL（Courant-Friedrichs-Lewy）≤1可以在生成的单元格网格内捕获流动特性。本研究将时间步长$ \Delta t $设置为0.02 s，采用一阶时间离散格式进行收敛，该收敛方式在多个研究证实是可行的^[18]。

3 数值方法验证 3.1 物模实验对比验证

为了对本文使用的CFD方法所计算出的结果进行可行性评估，本文对DTC箱体建模后满载情况下开阔水域的h/T为1.2、2.0两种工况的结果与弗兰德斯水力学研究所（FHR）^[19] 的物模实验结果进行对比验证。如图4所示，数值模拟结果与物模实验数据吻合良好，不仅验证了本文数值模拟计算方法的可行性，也验证了箱体建模后的DTC船模可用性。当水深吃水比h/T=1.2时，重心下沉量S_COG=15.60 mm，与FHR物理模型实验结果的偏差为5.45%；在h/T=2.0工况下，重心下沉量S_COG=6.32 mm，与实验结果的误差为1.86%，S_AP=5.4 mm，与FHR的物模实验结果误差仅为5.59%，S_FP=8.336 mm，与FHR的物模实验结果误差仅为0.7%，这2种工况下实验误差均小于10%，在可接受的误差范围内。

3.2 网格收敛性分析

为评估网格尺寸与时间步长引起的空间离散误差和时间离散误差，本文对吃水14.5 m、h/t=2.0的工况进行网格收敛性分析。使用网格细化比$r= \sqrt{2} $生成3组网格，3组网格采用同一CFL条件，且严格保证$ CFL< 1 $，其表达式见式(6)。模拟结果见表2，可知当网格数为1.5×10⁶时结果已较为精确，为节省算力资源，本文所有工况均采取此网格结构。

$ CFL=\frac{v\Delta t}{\Delta x}。$

(6)

式中：$ v $为船速；$ \Delta x $为时间步长；$ \Delta t $为网格间距。

4 仿真结果与分析

铜鼓航道水深为17.5 m，以赤湾理论最低潮面起算，年最高潮位3.80 m，因此本文航道水域计算域工况的水深与吃水比h/T设置为1.2、1.3、1.4。根据铜鼓航道的通航管理规定，船舶在7级风以上和台风期不能通航，因此本文数值模拟的工况风速设置为6级风。根据该水域风况统计资料风向设为西风。铜鼓航道落急横向流速为0.22～0.56 m/s，涨急横向流速为0.12～0.24 m/s，流况设置2种，分别为东南流，流速为0.514 m/s；西流流速为0.617 m/s。基于此，本文可总体划分2种工况：1）水深吃水比1.2、1.3、1.4下的船受无风无流、受东南流、受东南流和西风的状态；2）水深吃水比1.2、1.3、1.4下的船受无风无流、受西流、受西流和西风的状态。

基于上述工况的描述，主要研究2种工况下船舶在铜鼓航道通航时在风流作用、流作用、无风无流3种状态下所受横向力、偏航力矩和下沉量等3种船舶运动姿态。因此，每种工况具有3种风流状态，为表述方便，下文统一状态1为无风无流状态，状态2为无风有流状态，状态3为有风有流状态。

4.1 船体横向力对比分析

如图5所示，船在3种状态下均会产生横向力和偏航力矩，同时对比2种工况中的状态2、状态3与状态1发现，受到风的影响时，横向力增加较为缓慢，因此流对船舶偏转作用特别显著，风对船舶的偏转影响次之。在无风无流状态下，船舶的横向力会随着水深的变浅而逐渐增大，即随水深吃水比（h/T）的减小呈递增趋势，在有风有流的状态下，船舶所受横向力随着水深的减小呈非线性成倍增加的趋势。无风无流的状态下，船舶在h/T=1.4下所受横向力为−0.167 N，随着水深的减小，船舶在h/T=1.3，1.2下产生的横向水动力分别增大了78.4%，151.5%，由此可见，当船舶处于浅水状态（1.2<h/T<1.5）时，流对船舶的操纵性将受到较明显的影响，船舶横向水动力会随着水深的减小而显著增大。当流自左舷正横前来，风自左舷正横后来（工况1）时，风流造成的横向力会形成“岸推效应”，使船舶偏离航道中心线航行，可以发现，船舶在h/T=1.4下所受横向力为−5.78 N，船舶在h/T为1.3、1.2下产生的横向水动力分别增大了41.2%、137.2%，可见，当受斜流和斜风作用时，由水深造成的横向水动力成倍增大，船舶的操纵性受到明显的影响，航向稳定性受到影响，对船舶保向、船舶靠离泊、船舶回转掉头的困难增加。当流自右舷正横前来，风自左舷正横后来（工况2）时，对比工况1的流速增加了0.103 m/s，在h/T为1.2、1.3、1.4下流产生的横向水动力对比工况1下横向水动力11.775、6.476、4.503 N分别增加了26.11%、27.51%、26.07%。可见，流速对船舶保向性有较大的影响，流速越大，船舶的横向漂移速度也越大，船舶在航行过程中更容易偏离规定的航线。对比2种工况下船舶所受风流状态的横向力对比发现，风对船舶造成的横向力不会随着水深的减小而增大，因此风对船舶的横向力主要取决于风速、风向，在风的作用下，水深变化对船舶保向性影响较为有限。

4.2 船舶所受偏航力矩分析

根据图5和图6可知，船舶在工况1下，船舶在受风流耦合作用下的横向力大于船舶在流作用产生的横向力，但所受偏航力矩的大小情况却与之相反，船舶只受流作用产生的偏航力矩大于在风流耦合作用下的偏航力矩；在工况2下，船舶受横向力、偏航力矩的大小情况与工况1相反，在风流作用下横向力虽然增大但是偏航力矩却减小了。

2种工况下船舶所受风的大小和方向一致，故2种工况下风对船的作用一起分析。船舶在风中航行时，当风自左舷正横后来时，风力作用点靠近船尾，风力会在船舶的左舷正横方向上产生一个向右的横向分力，这个横向分力会试图将船舶的左舷向右推，在这个横向力的作用下，船舶对水产生一个相对运动，从而产生一个使船舶向左偏航的水动力矩，船首向上风偏转。

船舶在航行时受到斜流的作用，会使船舶受到横向力的作用，从而产生偏航力矩。在工况1下，流自左舷正横前来，水流会在船舶的左舷正横方向上产生一个向右的横向分力。这个横向分力会试图将船舶的左舷向右推，从而产生一个使船舶向右偏航的力矩。在工况2下，当流自右舷正横前来时，水流会在船舶的右舷正横方向上作用，水流的作用点主要集中在船舶的右舷侧，因此水流产生的横向力作用在船舶的右舷侧，这个力会使船舶的尾部向左偏航，从而产生一个使船舶向左偏航的力矩。

综上，在工况1下，尽管风流对船舶造成的横向力一致，船舶的偏航力矩是流和风产生的力矩的矢量和。由于流和风产生的力矩方向相同，它们会相互叠加，最终的偏航力矩会增大；在工况2下，流和风产生的偏航力矩方向相反，流产生的偏航力矩使船舶向左偏航，而风产生的偏航力矩使船舶向右偏航，因此这两个力矩会相互抵消一部分。由于流产生的偏航力矩大于风产生的偏航力矩，船舶最终向右偏航。

4.3 船体下沉量对比分析

船舶在航行过程中，船首部分将迎面水流沿舷侧向船尾方向分流，并在尾部区域形成汇聚流场。这种流动特性导致船体首尾两端压力升高，而舯部区域压力降低，形成纵向压力梯度，同时在浅水中，水流主要被限制于水平方向流动，导致船体下方动压升高，直接增大下沉量，且随着船速的增大，下沉量也会增大。根据图7可知，本文数值模拟中船顶流航行，船速“相对”减小，受流的影响导致“船速降低”的下沉量却大于船在浅水中的下沉量，这是由于斜向来流冲击船首，使船舶两侧的流速分布不对称，进一步产生横向压力差，该压差导致船体向下流侧倾斜，进而增大局部吃水。同时斜流在船首处产生局部高压区，阻碍水流顺畅通过舷侧，使船底流速进一步增加，舯部低压区扩展，下沉量增大。这也是船受斜顶流的作用导致船体的下沉量增大的原因。另外可以发现随着水深吃水比（h/T）的减小也会放大斜流对下沉量的影响，同时船在2种工况下航行，浅水效应引起的下沉量通常更大，斜顶流会额外增加下沉量，但其贡献量级低于浅水效应。值得注意的是，风力也会在船舶的纵向方向上产生一定的力，但是对船体的下沉量的影响较小，更多是体现在船体的横向力和偏航力矩上。

4.4 船体底部压力分析

为更清晰地理解风对船舶的下沉量的影响，本文选取工况1（流自左舷正横前来，风自左舷正横后来）下h/T=1.2 时DTC受流、风流作用下在铜鼓航道航行时进行船底部动压力分析，图8为去除船底静压后的船底部动压力分布情况。可知，这2种状态均在船首、船尾处捕捉到了动压力，反应了浅水效应和斜顶流作用下船首尾的压力梯度变化，同时在船首处捕捉到了更大的动压力，这也解释了斜流在船首处产生局部高压区，会阻碍水流顺畅通过舷侧，使船底流速进一步增加，船体下沉这一现象。由于图8中顶部来流的状态一样，因此船首底部的最大动压力相等；图8(b)中船尾处来风，但是对船尾底部压力的变化几乎无影响，这也解释了风对船体下沉量和纵倾影响较小的原因。

5 结　语

1）流速增大导致船舶横向力线性增加直接加剧船舶向航道一侧的横向漂移，这种漂移可能导致船舶与岸壁距离缩小，引发岸壁效应，增加碰撞岸壁的风险。流速方向会改变偏航力矩的方向：左舷来流产生向右偏航力矩，右舷来流产生向左偏航力矩。同时，流速增大使偏航力矩幅值显著变化，导致船舶航向难以稳定，需频繁调整舵角补偿力矩。

2）数值结果表明，风对船体的下沉量和纵倾影响较小。风恒定时，其对船舶的横向力和偏航力矩不会随着水深的减小而逐渐增大，受水深影响较小。

3）当流速增加0.103 m/s时，在水深吃水比（h/T）为1.2、1.3和1.4的工况下，横向水动力相较于工况1的11.775、6.476、4.503 N分别提升了26.11%、27.51%和26.07%。结果表明，流速对船舶保向性具有显著影响，其变化趋势表现为流速增大时，横向漂移速度随之增加，导致船舶在航行过程中更易偏离预定航线。

4）流速对横向漂移和偏航力矩的具体影响幅度（如流速每增加 0.1 m/s，横向力增幅约26%），可为大型集装箱船在不同流速的航行环境下提前调整航速（如降低航速以减小水动力）、优化舵角（如预判偏航方向并提前施加反向舵力），避免因保向失效导致危险。

此外，本文所用船模按照 1∶89.11的比例进行建模，但未考虑尺度效应，后续可进一步深入探究这一影响，研究大型集装箱船在铜鼓航道时通航需要多大的舵角可以克服风流的影响，未来还可以考虑在不同的船速下，风流造成的横向力和偏航力矩及岸壁效应的影响。为未来更多的大型集装箱船能够在铜鼓航道安全通航提供参考，保证航行安全。