0 引　言

通常船舶螺旋桨的设计只关注直航情况下的推进载荷，如推力、扭矩。然而，船舶在回转或遭遇横流、横风时，螺旋桨即在斜流中工作。斜流不仅影响螺旋桨的推进性能，而且还会产生较大侧向力和弯矩，改变船尾轴承载荷。当漂角较大时，侧向载荷可能超过设计值，从而导致轴系统失效。Amini等^[1]在实验中表明，斜流中螺旋桨侧向载荷引起的轴承载荷可达到只有螺旋桨重量引起的载荷的3倍，这容易造成船舶尾轴系统的过载损坏，影响船舶的安全性。近年来，随着船舶设计越来越轻量化，船体结构刚度下降，尾轴系统在侧向载荷作用下易出现大幅变形的情况，甚至出现轴承损坏的事故。因此，斜流中螺旋桨侧向载荷的预报显得更加重要。

目前，学术界已对敞水中螺旋桨斜流的水动力学性能开展了一些研究，如Dubbioso等^[2]，Wang等^[3]，Hou等^[4]。研究表明斜流中螺旋桨载荷会增加，其大小与漂角的变化有关，漂角较大时，会产生较大的侧向载荷，因此在轴系统设计中重视螺旋桨的侧向载荷十分必要。然而，实际上螺旋桨工作在复杂的船体尾流中，这使得其水动力载荷与敞水情况差异明显。根据Xing等^[5]的研究，船体尾流受斜流的影响较大，由斜流引起的船后螺旋桨侧向载荷的变化将更加复杂。因此，近年来船后螺旋桨的斜流效应得到了更广泛关注。Muscari等^[6]数值分析了斜流中双螺旋桨船尾部附体对侧向载荷的影响，详细讨论了螺旋桨侧向载荷下，螺旋桨-船体伴流的相互作用。Zhang等^[7]基于商业软件计算了斜流条件下潜艇模型的螺旋桨水动力。计算结果表明，随着攻角的增大，螺旋桨的推力和力矩先减小再增大，为潜艇机动性研究中螺旋桨的水动力预测提供了一种方法。Ortolani等^[8]的研究表明，螺旋桨的侧向载荷比轴向载荷对空泡更敏感，而实尺度下的空泡会导致侧向载荷的大量增加，基于半经验方法对螺旋桨载荷的尺度效应进行了初步分析。Ortolani等^[9]进行了双螺桨船纯斜流试验，对纯斜流和自由转弯的结果之间的差异进行了讨论。结果发现，由于迎风螺旋桨和背风螺旋桨的流入特性不同，其差异有所不同。Yao等^[10]在开源平台OpenFOAM上，通过RANS模拟研究了螺旋桨在斜流中的水动力性能。采用滑动网格的方法计算了螺旋桨的旋转运动。在无空化或弱空化条件下，实验结果与实验数据基本一致。Sun等^[11]对螺旋桨旋舵系统的KCS船模进行了数值预测，结果发现，随着漂角的增大，阻力、侧力和偏航力矩的系数不断增大，且漂角越大，螺旋桨盘的整体均匀性越差，漂角的螺旋桨推力和转矩波动系数峰值大于直线导航，负漂角大于正漂角。相对而言，关于斜流工况下螺旋桨特性的研究较少，未形成统一的预报方法，尤其是实尺度下，螺旋桨斜流载荷的预报还没有有效的换算方法。

螺旋桨侧向载荷估算是船舶尾管轴承载荷计算和斜流轴系统结构安全性能分析的基础。由于缺少这方面的研究数据，工程设计中根据经验简单估算，增加了实船运营中轴系出现安全事故的风险^[12]。

随着CFD的发展，实尺度计算技术逐渐成熟，为一些复杂流动问题的尺度效应研究提供了解决途径。Wang等^[13]研究了一艘集装箱船在直航时的实尺度载荷和激振力，Zhang等^[14]分析了一艘散货船在斜流中，螺旋桨侧向载荷的尺度效应。本文基于CFD数值模拟针对集装箱国际标模KCS开展实尺度船舶斜流中的水动力计算研究，重点分析实尺度下，斜流对螺旋桨载荷的影响，为船舶螺旋桨设计和轴系安全计算提供参考。

1 几何模型与网格 1.1 KCS 介绍

KCS采用单桨推进，实长230 m，设计航速为24 kn，几何形状如图1所示。螺旋桨为5叶桨，主要尺度如表1所示。作为一艘国际标模，KCS已被很多学者研究过，如阻力计算^[15]、自航模拟^[16]，波浪中运动^[17-19]、操纵运动^[20]等。

1.2 计算域及网格

本文数值模拟采用商业CFD软件STAR-CCM+，网格基于软件中裁剪单元方法生成，该方法能够产生非结构化的六面体网格。图2为实尺度船的计算域网格划分。总网格数为835万，船体网格数为482万，螺旋桨网格数为353万。为捕捉更精准的计算结果，分别在船体自由表面、船尾、螺旋桨及其导、随边等部位对网格加密。螺旋桨旋转采用了滑移网格方法。

2 数值方法

流动模拟采用非定常的雷诺平均Navier-Stokes(URANS)方法，控制方程可用不可压缩流体的张量表示如下：

$ \frac{{\partial\bar{u}}_{\rm{i}}}{\partial{{x}}_{{i}}}=0，$

(1)

$\begin{split} \frac{\partial }{\partial t}\left(\rho {\bar{{u}}}_{{i}}\right)+&\frac{\partial }{\partial {x}_{j}}\left(\rho {\bar{u}}_{{i}}{\bar{u}}_{{j}}\right)=-\frac{\partial \bar{P}}{\partial {x}_{i}}+\\ &\frac{\partial }{\partial {x}_{j}}\left[\mu \left(\frac{\partial {\bar{u}}_{i}}{\partial {x}_{j}}+\frac{\partial {\bar{u}}_{j}}{\partial {x}_{i}}\right)\right]+\frac{\partial }{\partial {x}_{j}}\left(-\rho \overline{{{u}}_{{i}}'{{u}}_{j}'}\right){+\bar{f}}_{i} 。\end{split}$

(2)

其中： $ i $ 和 $ j $ 为坐标分量， $ u $ 为速度， $ \rho $ 为流体密度， $ \mu $ 为流体的粘度， $ P $ 为压力， $ f $ 为体积力， $ \rho \overline{{u}_{i}'{u}_{j}'} $ 为雷诺应力项。

为了封闭控制方程，采用了湍流模型SST k-ω。SST k-ω模型将k-ε和k-ω模型结合在一起，使其非常适合于模拟具有分离和强逆压梯度的流动。k和ω的输运方程如下式：

$\begin{split} \frac{\partial }{\partial t}\left(\rho k\right)+& \nabla \cdot \left(\rho \bar{u}k\right)=\nabla \cdot \left[\left(\mu +{\sigma }_{k}{\mu }_{t}\right)\nabla k\right]+\\ &{G}_{k}-\rho {\beta }^{*}{f}_{\beta *}\left(\omega k-{\omega }_{0}{k}_{0}\right)，\end{split} $

(3)

$ \begin{split}\frac{\partial }{\partial t}\left(\rho \omega \right)+& \nabla \cdot \left(\rho \bar{u}\omega \right)=\nabla \cdot \left[\left(\mu +{\sigma }_{\omega }{\mu }_{t}\right)\nabla \omega \right]+\\ &{G}_{\omega }-\rho \beta {f}_{\beta }\left({\omega }^{2}-{{\omega }_{0}}^{2}\right) 。\end{split}$

(4)

其中， $ {{\sigma }}_{{k}} $ , $ {{\sigma }}_{{\omega }} $ 为模型系数。

$ {\sigma}_{{k}}=0.85{{F}}_1+1-F_1\text{，}{\sigma}_{\omega }=0.5F_1+0.856(1-F_1)。$

(5)

式中， $ {{F}}_{{1}} $ 为混合因子，计算中根据局部流动情况来确定具体值。

式(3)和式(4)中的 $ {G}_{k} $ 和 $ {G}_{\omega } $ 为湍流产生项。 $ {f}_{\beta *} $ 为自由剪切流引起的修正， $ {f}_{\beta } $ 为涡伸缩引起的修正。 $ {k}_{0} $ 和 $ {\mathrm{\omega }}_{0} $ 为周围环境湍流值。以上系数计算如下：

$ {G}_{k}={\mu }_{t}{f}_{c}{S}^{2}-\frac{2}{3}\rho k\nabla \cdot u-\frac{2}{3}{\mu }_{t}{(\nabla \cdot u)}^{2}，$

(6)

$ {G}_{\omega }=\rho \gamma \left[\left({S}^{2}-\frac{2}{3}{(\nabla \cdot u)}^{2}\right)-\frac{2}{3}\omega \nabla \cdot u\right] 。$

(7)

其中： $\nabla \cdot {u}$ 为速度散度， ${S}$ 为应变率张量的模。

$ {f}_{\beta *}=\left\{\begin{array}{l}1，\qquad\qquad\quad\;\;\;{\gamma }_{k}\leqslant 0\\ \dfrac{1+680{{\gamma }_{k}}^{2}}{1+400{{\gamma }_{k}}^{2}},\qquad\;{\gamma }_{k}\geqslant 0\end{array} \right.$

(8)

其中：

$ {\gamma }_{k}=\dfrac{\nabla k\cdot \nabla \omega }{{\omega }^{3}} ，{{f}}_{\mathrm{\beta }}=\frac{1+70{\mathrm{\gamma }}_{\mathrm{\omega }}}{1+80{\mathrm{\gamma }}_{\mathrm{\omega }}} _{ }，$

(9)

式中， $ {{\gamma }}_{{\omega }}=\dfrac{\left({W}\cdot {W}\right):{S}}{{\left({{\beta }}^{{*}}{\omega }\right)}^{3}} $ 。其中： $ {W} $ 为速度的旋转张量； $ {\ \beta }^{*} $ 为系数，大小由流动状态确定。

压力-速度耦合采用了SIMPLEC算法，时间离散采用二阶隐式格式，对每个时间步长设置5次内部迭代。对流项采用了二阶精度迎风格式，根据VOF法捕获自由表面。

VOF 控制方程为：

$ \frac{\partial }{\partial t}\left(f\right)+\frac{\partial }{\partial {x}_{i}}\left(f{u}_{i}\right)=0 。$

(10)

其中，f为流体体积分数。

根据ITTC^[21]的CFD计算指导意见，螺旋桨旋转运动问题的时间步长应至少在一个螺旋桨旋转周期内有200个时间步，本文计算中设置时间步长为： $ \Delta t= 1/\left(200n\right) $ ，n为螺旋桨转速，满足ITTC指导意见里的要求。为防止船体兴波在计算域边界处出现反射之类的非物理数值问题，在压力出口附近设置了消波功能，消波区长度为0.5L_PP，消波方法为STAR-CCM+自带的阻尼消波法。

3 计算可靠性分析

由于KCS还没有建造成实船，因此目前还没有KCS的实船数据。以其他实尺度计算文献为参照，分析本文的计算可靠性。计算工况选择设计航速24 kn时的自航工况，计算结果如表2所示。其中，EFD数据是根据模型试验外插得到^[22]。从计算结果看，本文得到的螺旋桨转速、推力和扭矩整体上与试验预报结果和其他文献结果比较接近，吻合较好。

4 计算结果分析 4.1 螺旋桨载荷

漂角的定义如图3所示。横向水流来自右舷侧为正，其中 $ \ \beta $ 为漂角，U₀为横向水流， $ {U}_{X} $ 为船长方向速度，U_y为船宽方向速度，分别对漂角为 $ \ \beta $ =0°，±5°，±10°，±15°，±20°的情况展开计算，获取桨在X、Y、Z这3个方向的受力、力矩。图4给出了坐标轴位置和方向，取桨盘中心处为原点，船尾指向船首为X轴正向，右舷指向左舷为Y轴正向，船底指向甲板为Z轴正向。

图5为不同漂角 $ \ \beta $ 下螺旋桨的受力。由图5（a）和图5（b）可知， $ \ \beta $ 为正时，推力和扭矩大小略有减少， $ \ \beta $ 为负时，推力和扭矩大小略有增加。原因主要是当船体阻塞横向水流后，存在一个绕着船体底部的流动，使桨盘面下半面的横向速度变化较大，对于正向漂角而言，下半桨盘面的绕流起到减小桨叶载荷的作用，而负向漂角则相反。由图5（c）和图5（d）可知，螺旋桨横向、垂向的受力因漂角 $ \ \beta $ 不同变化较大， $ \ \beta $ 越大，横向受力明显增大，会影响船体的操纵运动。此外， $ \ \beta $ 为−20°时，垂向受力明显增大，其原因是横向水流使船体伴流产生偏移，使得桨盘处的伴流场左右严重不对称。

图6为横向和垂向载荷与轴向载荷比率随漂角的变化。横向与轴向载荷比率随漂角成非线性变化，漂角角度越大，横向载荷占比越大，对应20°漂角，占比达68.72%；漂角为−20°时，垂向与轴向载荷比值达72.91%。

4.2 桨叶表面动压力

各漂角下螺旋桨表面动压力分布如图7所示。可看出，随着漂角的增加，桨盘面下半盘面载荷逐渐减小，而上半盘面的载荷逐渐增加，这是由于随着攻角的增加桨叶来流攻角在下半盘面逐渐减小而在上半盘面逐渐增加。这一规律能够比较清晰地从单独一个桨叶的推力载荷中看到，如图8所示，其中桨叶的相位角定义如图9。可看出，在同一相位角处，攻角的增加显著改变了桨叶载荷，尤其是在0°相位附近。此外，图8也显示漂角的存在很大程度地改变了桨叶载荷的变化规律，对比−20°和20°工况下的桨叶载荷随相位角变化曲线可看到二者存在较大的差别。这种差别对螺旋桨的设计有较大意义，因为可能会造成螺旋桨空泡形态的较大变化，进而影响剥蚀、激振力等现象。从桨叶相位角看，不同相位处的桨叶载荷对漂角的敏感性不同，相对而言，桨盘面最低位置附近（0°相位处）对漂角更加敏感，而桨盘面最高位置附近（180°相位处）对漂角敏感度较低。这主要是因为桨盘面下端附近船体遮蔽效应较弱，横向水流影响更加显著。从图10的流线也可看出，水流绕过船体尾轴底部时存在偏转，这种流线的偏转在螺旋桨轴下方流动区域较上方流动区域更加明显，当这种绕流的流动方向与螺旋桨旋转方向一致时，减低螺旋桨桨叶载荷，反之则增加螺旋桨桨叶载荷。

4.3 尺度效应分析

对于船舶流动问题，实船尺度和模型尺度的雷诺数相差较大，这导致船体和螺旋桨受力在模型尺度和实尺度下不一致，即使以无量纲的形式表示仍然存在较大差异。对于快速性问题，国际拖曳水池（ITTC）给出了建议的换算方法，在工程界应用较广。一些学者曾通过数值模拟分析尺度对船体尾迹和螺旋桨轴向负载的影响，如Castro等^[22]，Wang等^[23]。这些研究也表明由于雷诺数的变化，模型与实尺度结果存在一定差异。

相对于传统的快速性问题，斜流下的螺旋桨受力则尚未有成熟的换算方法。本文采用的数值计算方法与文献[24]采用的方法相似，两者有可比性。图11为在不同尺度和漂移角度下的螺旋桨载荷。结果表明，实尺度下桨的推力对漂角更敏感，横向载荷变化规律对漂角不敏感。相比之下，垂直载荷对尺度较为敏感，且尺度影响的大小与漂角有关，漂角增大时，实尺度与模型尺度结果差异变大。

5 结　语

本文基于CFD数值模拟国际标模KCS实尺度下斜流工况的流动问题，分析斜流对螺旋桨载荷、桨叶表面压力分布、船体周围流动的影响。计算采用RANS方法，结合SST k-ω两方程湍流模型和滑动网格方法，在0°漂角下通过对比文献结果验证了计算的可靠性。为分析斜流的影响，计算考虑了±20°，±15°，±10°，±5°以及0°漂角。结果表明，斜流对螺旋桨水动力载荷影响较大，斜流对推进力及力矩的影响相对漂角不对称，主要是因为来自左右舷侧的横向水流对螺旋桨盘面载荷的影响不同。斜流的存在导致螺旋桨受到较大横向力和力矩，方向和大小大体上与漂角呈反对称分布。此外，通过单独桨叶的推力变化曲线分析了漂角对不同相位角度处螺旋桨载荷的影响规律。结果表明，在桨盘面下端附近桨叶载荷更为敏感于漂角变化。此外，计算尺度对斜流下的螺旋桨载荷系数有一定程度影响，相对而言轴向载荷和垂向载荷系数的尺度效应较大。