在碎冰区航行的船舶，其螺旋桨会受到水下碎冰的影响，发生阻塞和碰撞等形式的干扰，导致螺旋桨所受载荷将远大于敞水工况^[1], 从而导致螺旋桨强度受到较大的影响^[2]。对于螺旋桨所受冰载荷的作用，可以分成接触载荷和非接触载荷^[3]，众多学者对此进行了大量的数值计算与试验研究。对于非接触载荷，Liu等^[4]使用面元法代码PROPELLA分别在敞水状态下和阻塞流中进行了旋转的螺旋桨的数值模拟，计算了不同的冰桨间距下螺旋桨表面受力情况。常欣等^[5]使用Fluid计算了在冰桨间距变化过程中的螺旋桨水动力性能。文献[6]在空泡水桶中开展冰阻塞模拟实验，研究冰阻塞下冰厚比、冰桨间距及空泡数对螺旋桨性能的影响。在接触载荷上，Ye等^[7]等根据近场动力学理论，采用自编程方法计算了在冰桨碰撞铣削过程的螺旋桨桨叶的载荷情况；王建强等^[8]利用CCS和IACS规范, 对冰载荷下螺旋桨叶片进行强度分析, 计算水动力载荷和冰载荷对螺旋桨强度的影响。以上的这些研究多关注于螺旋桨载荷，但是在冰桨干扰过程中，影响是相互的，冰的运动轨迹会影响冰桨碰撞概率，碰撞过程桨叶和桨轴以及主机都将受到较大的影响^[9]，并且碰撞会加大桨的振动与噪声，综合来看，研究冰的运动轨迹对于冰桨碰撞至关重要。Veitch^[10-11]对冰的受力情况展开研究，总结出冰的运动受冰形状、重力、浮力以及附加质量力等的影响；并且约束力主要来自于冰块之间，冰船以及冰桨之间的接触力。为了研究方便，Veitch简化冰桨干扰模型，分析切削过程中冰的受力。在冰桨干扰过程中冰同样会受到2种力：一种是接触力，这是冰桨碰撞所引起的，并且对于冰和桨来说接触力是相互的；另一种是非接触力，主要是来自于螺旋桨抽吸作用下的水动力，这也是本文的主要研究内容。本文建立了特定的冰桨干扰模型。采用重叠网格方法对流域生成网格进行计算，在此基础上，变化螺旋桨转速和冰桨径向间距来进行变工况计算。

1 运动模型的建立 1.1 基本模型的设置

螺旋桨计算模型选用Icepropeller Ⅰ桨^[12]，螺旋桨直径D=0.2 m, 叶数为4叶，螺距比为0.78，毂径比为0.3，盘面比为0.67，纵倾角为10°。如图 1所示，将坐标原点选在螺旋桨的中心处，沿桨毂长度方向设为X轴方向，桨叶的母线设为Y轴方向，X、Y、Z轴遵从笛卡尔坐标系。由于桨盘面是在YOZ平面内，后文中冰块沿Y轴方向的运动称为水平运动，沿Z轴的运动则称为竖直运动。模型冰设置为边长0.1D的正方体形状，其中D为螺旋桨的直径。冰与螺旋桨的相互作用模型如图 1所示。

1.2 网格的划分

本文主要研究外流场问题，所以在建立外流域时也可认为外流域对内部流场没有影响。即假定外流域无限大，具体计算区域如图 2所示，其中外流域长12D，直径为6D，入口设为速度入口，出口设为压力出口。

进行网格划分前需把各区域的边界条件进行设定，其中在冰和螺旋桨的外部区域各新建一个加密区，并设置为重叠网格，螺旋桨表面、冰表面和外流域旋转表面均设置为壁面边界。设置完成后选取网格类型包括切割体网格和棱柱层网格，为了减少计算量、节省计算时间、提高服务器的应用效率，把对计算域没有影响的外流域壁面网格中的最小相对尺寸及相对目标尺寸都设为120%，网格划分如图 3所示。

为了得到更加精细的流场，需对桨叶及周围进行单独加密处理，并且在冰桨之间的运动区域也要进行网格加密，方便重叠网格间的数据传递，网格大小数量级需与冰桨重叠网格相同，具体加密情况如图 4所示。

1.3 物理模型的选取

本文计算的主要目的是得到桨前冰的运动情况，在选择物理模型时要尽量考虑到真实工况下的冰桨相互作用关系，因此需进行隐式非定长计算，选择的物理模型包括液体、恒密度、湍流和k-ω模型，其中k-ω模型能够较好地处理不同压力梯度下的边界层问题。

计算过程中的冰为自由运动，在6个自由度上均无约束，需利用动态流体相互作用(dynamic fluid body interaction, DFBI)模型计算。DFBI模型的原理是通过计算刚体在流体中的受力以及刚体与刚体之间的相互作用力推导其运动轨迹，且DFBI模型的影响因素有很多，包括物体的质量、物体所受到的惯性力矩、通过积分流体压力而产生的力、自定义的力和力矩，DFBI运动的求解器是在全部6个自由度下进行计算的。在STAR-CCM+中DFBI属性是模拟树中一个独立的输入窗口，在模拟计算时至少一个区域要被指派到DFBI对象中。求解器设置过程中，零部件选为冰，质量为0.006 4 kg，质心坐标(-0.19 m，0，0.09 m)，对于X轴的惯性距(5.227×10^-5 m，5.227×10^-5 m，5.227×10^-5 m)，初始速度为1.6 m/s。

2 数值计算结果分析

选取冰在桨前的轴向位置D_A=0.19 m，且保持不变，通过改变径向位置D_R和进速系数J来变更工况，并对比不同工况下的计算结果总结出相应的规律。计算时选取了3种进速工况，5种径向位置。3种工况下的进速系数J分别取0.2、0.4、0.6，5种径向位置分别为0.4R、0.5R、0.8R、1.0R和1.2R，来流速度以及冰的初始速度V_A=1.6 m/s。

2.1 干扰工况下冰的受力分析

改变冰运动状态的力主要是由螺旋桨的抽吸作用所产生的。螺旋桨在旋转时对周围的流体产生一定的扰动，且越靠近桨叶部分周围流体的速度越大，这就造成了冰前与冰后的流体速度不一致，进而改变了冰的轴向速度。而螺旋桨的抽吸作用还会引起周围流体的转动，因此对冰在周向上的运动状态也将产生一定的影响。

通过逐一计算，获取了从桨前同一轴向距离开始运动时冰的受力情况，轴向、周向和径向的受力具体计算结果如图 5~7所示。

由图 5可知，进速系数越小，冰块沿X轴正方向上的受力的峰值越大，且都出现在桨盘面处。这是由于螺旋桨的抽吸作用使得桨盘面处的流体流速增大，引起了冰块前后表面的速度差达到最大，因此对冰块的作用力也达到最大值。当冰块运动到桨后方较远处时(桨后轴向距离约4R)，冰在轴向上的受力趋于零，这是由于桨的抽吸作用衰退，使得流体的流速趋于初始值。

进速系数J=0.2时，受力的最大值点出现在D_R=0.5R处；J=0.4时，最大值点出现在D_R=0.8R处；J=0.6时，最大值点出现在D_R=R处。因此可以得出随着进速系数的增大，5种工况下，冰块受力的最大值点出现在径向距离较大处。而同一进速系数下，径向距离越大，受力的峰值点出现的越早；同一径向距离时，则进速系数越大，受力的峰值点出现的时间越晚。

在图 5中之所以会出现负值，是因为开始时螺旋桨的抽吸作用不明显，此时螺旋桨相当于一个圆盘面，对前方来流有一定的阻塞作用，也称为壁面效应，因此冰块受到水流的反作用力而出现负值。

图 5(b)、(c)中的D_R=0.04与D_R=0.05出现截止点是因为此时冰桨发生了碰撞，计算停止；而图 5(a)没有这种现象，是因为进速系数J小，螺旋桨的转速大，冰受到的作用力也较大，有被“甩”出去的趋势，因此冰桨发生碰撞的时间较晚。

图 6描述了冰桨相互干扰时冰块在Y方向上的受力情况，引起水平方向受力的原因是螺旋桨的旋转带动了周围水流转动，转动方向与螺旋桨旋转方向相同，而冰块在水流的冲刷作用下开始有水平方向上的速度。由图 6知，冰块在Y方向上的受力随着冰桨间距的减小而增大，且在桨盘面处达到最大值，在脱离桨盘面后又逐渐减小，这是因为越靠近桨盘面处，螺旋桨的转动作用使得水流获得的周向诱导速度越大，作用在冰块上的力也就越大。当水流过螺旋桨后，作用在流体上的外力等于零，此时冰不再受到向外的水平方向的力，而冰的速度大于水流的速度，因此受到阻力，即水平力出现负值。等到冰块脱离干扰区后，冰块不再受到抽吸作用，Y方向上的受力逐渐恢复到零。

同时，由图 6可以看出冰块在Y方向受力的最大值点，也随着进速系数的增大，出现在径向距离较大的工况下。图 6(b)、(c)中的D_R=0.04与D_R=0.05出现截止点的原因与X方向受力时出现截止点的原因相同，都是由于此时冰与螺旋桨发生了碰撞。

图 7反映了冰块沿竖直方向上的受力情况，受力的主要原因也是螺旋桨的旋转使周围流体产生了径向的诱导速度，使得冰块在旋转的同时产生了一定的离心力，具有脱离桨盘面运动的趋势。因此，在Z方向上产生了一定的受力，且力的大小与螺旋桨的转速有关，转速越大，产生的竖直方向力越大。由图 7可知冰块在Z方向的受力在桨盘面处达到最大值，同一进速系数下，不同径向距离的冰块在Z方向上的受力变化趋势大致相同；不同进速系数下，受力的最大值点出现在不同的竖直位置处，且进速系数越大，峰值点出现在的竖直位置越远。图 7中会出现正负值交替变换现象是因为冰块由于螺旋桨的旋转有被甩出去的趋势，紧接着由于水流的作用向下运动，即沿Z轴的负方向，因此受力出现负值。

2.2 干扰工况下冰的运动分析

根据之前计算得出的冰在各个方向的力，通过式(1)可以积分得出其位移变化：

$ S={{S}_{0}}+{{V}_{0}}t+\iint{f\left( t \right){{\text{d}}^{2}}t} $

(1)

式中:S₀表示冰在3个方向上的初始位移；V₀表示初始速度；f(t)表示上文所提到的冰的受力。

从图 8可以看出，多数情况下X方向的位移接近于一条直线。分析式(1)可知，在X方向上的初始位移S₀=-0.19 m, 初始速度V₀=0.16 m/s。可见，式(1)中的第3项积分项相对于V₀是一个小量，这从另一个侧面也反映了在冰块径向运动过程中螺旋桨抽吸作用影响不大，初速度占了主要的影响。

图 9反映了冰块在Y方向的位移变化，由图可知冰块在Y方向上的位移大致符合二次函数。引起Y方向上位移的原因是冰块受到切向力的作用，切向力以及流体阻力的合力给冰块一个沿切向的加速度，从而引起了冰块的位移曲线图接近二次函数。

初始时刻冰块在Y方向上位移的增加较缓慢，这是由于桨的抽吸作用在远离桨盘面处不明显；大约在0.1 s处位移迅速增加，此时冰桨间距小于R，螺旋桨的抽吸作用明显；而当冰块运动到桨盘面时，螺旋桨的抽吸作用达到最大，冰块在Y方向上的受力值也最大；在冰块越过桨盘面后，抽吸作用产生的力开始减小，冰块在Y方向上的位移变化速度也开始减小。图 9(b)、(c)中的截止点是由于此时冰桨发生了碰撞。

由图 10知冰块在Z方向上的位移变化与在Y方向位移变化基本相同，在初始时刻位移的增加较缓慢，大约在0.1 s处位移迅速增加，而当冰块运动到桨盘面时，螺旋桨的抽吸作用达到最大，使得冰块在周向上的受力增加，旋转过程中的离心力也增大；在冰块越过桨盘面后，抽吸作用产生的力开始减小，位移增量也开始减小，直到脱离干扰区。

图 10(c)中出现负值因为此工况下螺旋桨的转速较慢，对周围流体的扰动较小，流体带动冰块旋转时的角速度很小，冰块在短时间内相当于直线运动，离心力的作用几乎可以忽略，而桨的抽吸作用又使得远离桨盘面处干扰区域中的流体压强减小，给冰块一沿Z方向向下的作用力，冰块在逐渐靠近的过程中，离心力越来越明显，冰块的位移恢复到正常情况。

同一进速系数下，水平距离越小，冰块沿Z方向位移的变化越快，且达到的峰值越大；不同进速系数下，冰块在Z方向的位移大小，随着进速系数的增大而减小，且变化范围也随之减小，最后达到的稳定值也较小。

3 结论

1) 螺旋桨的抽吸作用会使冰块在各个方向上的受力改变量峰值出现在桨盘面附近。

2) 螺旋桨进速系数对冰桨碰撞概率影响较大，可以通过适当调节进速系数大小减小冰桨碰撞概率。

本文初步计算分析了桨前自由状态冰的受力和运动轨迹，得出了一些有意义的结论，但本文是针对模型桨开展的计算，是否与实桨工况下冰的受力和运动规律类似，接下来还需要进一步研究。