随轴系做复杂空间运动的船舶螺旋桨水动性能计算

邹冬林，张建波，田佳彬，塔 娜，饶柱石，4

（1.上海交通大学 振动、冲击、噪声研究所，上海200240；2.上海交通大学 机械系统与振动国家重点实验室，上海200240；3.武汉第二船舶设计研究所，武汉430205；4.高新船舶与深海开发装备协同创新中心，上海200240）

0 引 言

由于海洋来流的不均匀，船舶螺旋桨工作中产生静推力的同时也不可避免地产生脉动力。该脉动力传递到轴上引起轴系振动，进而引起船体振动并辐射噪声。过大的振动会导致轴系、船体等结构疲劳损坏，产生安全隐患。同时由振动而产生的辐射噪声也会影响乘员舒适性。因此在船舶推进轴系设计阶段，需要考虑如何减少螺旋桨脉动力的产生以及优化螺旋桨脉动力到轴系的传递途径等一系列问题。为解决这一问题，需要全面分析流体-螺旋桨-轴系间的耦合动力学特性。流体-螺旋桨-轴系耦合系统动力学分析本质上是一个流固耦合问题，且非常复杂。这是因为一方面不均匀流体引起的螺旋桨脉动力会导致轴系振动，而轴系振动又会带动螺旋桨振动，从而改变螺旋桨的脉动力。二者之间相互耦合，形成反馈机制。另一方面，即使螺旋桨工作在均匀来流中，由于螺旋桨或转轴不可避免地存在偏心，使得船舶推进轴系在绕自身中心线旋转的同时，又发生空间涡动（又称为进动），导致附在其上的螺旋桨也做复杂的空间旋转运动，从而使螺旋桨的进流速度每时每刻都在变化，进一步导致脉动力的产生。总之，螺旋桨脉动力和轴系振动是相互影响的，并且可能存在增强效应。根据相关文献资料，由于船体和轴系振动诱发了螺旋桨和伴流场耦合面的流体振荡，可使螺旋桨激励力产生的低频声辐射增加10～15 dB，因此非常有必要考虑轴系振动对螺旋桨水动性能的影响。目前针对流体-螺旋桨-轴系的双向流固耦合这一复杂问题的研究文献很少。国内外已有的文献均是把这一问题割裂为两个问题：流体-螺旋桨流固耦合动力学分析及流体载荷激励下的轴系振动。通常是先研究流体-螺旋桨流固耦合动力学特性，求解出螺旋桨对根部的激励力。研究中认为螺旋桨根部固定，忽略轴系振动的影响。然后把激励力作为集中载荷施加在轴系末端上（与螺旋桨连接的地方），进而研究轴系的振动响应等等，此时螺旋桨通常作为一个集中质量。该研究思路忽略了轴系振动对流场的反馈作用。当轴系振动幅值很小时，由此导致的螺旋桨水动力学性能变化可以忽略，这种分开处理的方法是满足工程要求的。而对于大型柔性船舶轴系，由于其细长比通常很小，轴系振动幅值不能忽略；同时由于大型船舶载重增加，通常选用大尺寸螺旋桨，由此导致的螺旋桨微小偏心都会产生很大的不平衡载荷，从而引起轴系剧烈振动。因此把流体、螺旋桨和轴系作为一个统一的整体去研究其动力学特性及其演化规律是必要的。

对这一复杂耦合系统的动力学分析首先要研究以下三个方面的子问题，如图1所示。第一个子问题是流体-螺旋桨双向流固耦合问题的研究；第二个子问题是螺旋桨脉动力引起轴系振动这一单向流固耦合问题的研究；第三个子问题是轴系振动导致的螺旋桨脉动力变化的研究。只有将这三个关键子问题摸透了，才能将所有子问题串联在一起，进而建立一个复杂的流体-螺旋桨-轴系双向流固耦合动力学模型。

第一个子问题主要研究螺旋桨叶片变形、桨叶动态应力及水动力预报等等。第二个子问题主要研究轴系振动响应、固有频率及轴承力等等。目前对第一个和第二个子问题均有大量研究文献[1-4]，而对第三个子问题的研究几乎没有。第三个子问题是第二个子问题的逆问题，其核心点是如何预报随轴系做复杂空间运动的螺旋桨的脉动力。而本文的研究重点集中在第三个子问题上，目的在于探索一种复杂空间运动的船舶螺旋桨水动性能计算方法。

图1 流体-螺旋桨-轴系的双向流固耦合研究的三个子问题Fig.1 Three key problems for fluid structure interaction of propeller-shaft sysytem

对于具有复杂运动的船舶螺旋桨水动性能预报的文献，国内外研究均很少。目前有少量文献研究具有简单运动形式的螺旋桨的水动性能。如陶尧森等[5]利用势流理论近似计算了波浪导致螺旋桨做升沉运动时的水动性能。喻欣[6]用Fluent软件计算了螺旋桨做升沉运动时的水动性能。黄苗苗等[7]也用Fluent软件做了类似的工作。郭春雨等[8]利用试验研究了具有升沉运动的螺旋桨的水动性能。Politis[9-10]利用不定常面元法预报了螺旋桨做升沉运动时的水动性能。Kinnas[11-12]等也利用不定常面元法预报了螺旋桨做升沉或者纵向或者扭转振动时的水动性能。在所有这些研究中，螺旋桨的运动形式均比较简单，基本都限定在单一的平动形式。而对于随轴系振动的螺旋桨而言，由于其有6个空间自由度，使得其运动轨迹复杂，不同自由度间还具有耦合现象。因此上面这些文献中的方法不具有通用性。

针对上述文献中的不足，本文构建了螺旋桨在流场中做复杂空间运动时的水动力预报数学模型，并利用不定常面元法求解了这一问题。利用本文的方法可以方便地研究螺旋桨随轴系振动而做复杂空间运动时的水动性能，分析轴系振动对螺旋桨水动力的影响，从而解决了上文提到的流体-螺旋桨-轴系流固耦合问题中的第三个关键子问题，为研究流体-螺旋桨-轴系双向流固耦合动力学特性做基础。同时本文的方法也可以用来预报船舶转弯、船舶升沉及纵摇振荡、螺旋桨启停及加速等复杂工况下的螺旋桨水动性能。

1 数学模型

目前已有多种成熟方法预报螺旋桨水动力，如升力线理论、升力面理论、面元法（又称边界元方法BEM）、CFD法等。升力线理论含有大量假设，因而只适用于轻载螺旋桨。升力面理论虽然可以较为准确地预报总推力和总扭矩，但由于没有考虑螺旋桨叶片的厚度，因而其预报的桨叶面压力分布不够准确。CFD方法尽管能准确计算桨叶表面压力分布，但是由于其计算量大，耗费时间长，因而应用并不广泛。相对升力面理论来说，面元法考虑了叶片厚度的影响，因此不仅可以准确计算总推力与总扭矩，而且可以较准确地计算桨叶表面的压力分布；相对于CFD方法来说，面元法计算量小，耗时短，同时相比于早期来说，面元法有了很大的改进，应用范围越来越广，比如可以考虑空化现象，结合边界层理论还可以考虑流动分离等复杂流动现象。因此本文采用基于扰动速度势的面元法来预报随轴系做复杂空间运动的螺旋桨的水动性能。

要建立随轴系振动的螺旋桨水动力预报数学模型，需要解决三个问题：一是建立轴系复杂空间运动的数学模型；二是建立具有复杂空间运动固体的流体边界条件；三是建立时变的尾涡数学模型。

1.1 轴系复杂空间运动的数学模型

图2 轴系振动与三个坐标系示意图Fig.2 The schematic of shaft vibration and three coordinate systems

坐标系oxyz与坐标系ox1y1z1的关系可以用投影角法表示[14]，如图3所示。设螺旋桨的自旋轴ox在x1y1和x1z1平面上的投影线与ox1轴的夹角为θy和θz。θy和θz称之为投影角。坐标oxyz与坐标系ox1y1z1的关系可表示为：

图3 投影角示意图Fig.3 The schematic of projected angle

式中：φ=ωt+θx，ω为螺旋桨的自转角速度。

从而坐标系OXYZ与坐标系oxyz的关系可表示为：

1.2 流体边界数学模型

式中：Vin表示螺旋桨的进流速度。当Vs与ωs为零［T］为单位阵时，其与传统的不考虑轴系振动的螺旋桨进流速度一致。r表示螺旋桨表面某点的位置向量，在oxyz中表示。上述式子假设螺旋桨叶片刚性，不发生弹性变形，如果叶片有弹性变形，则还需要考虑其变形速度。

假设：（1）流体无粘、无旋且不可压缩；（2）螺旋桨浸水足够深，即不考虑自由液面影响，同时流体域延伸到无限远；（3）不考虑空化影响；（4）轴系振动为小幅振动，叶片上没有涡的分离等复杂流动现象。取一足够大的外部控制面将其封闭在内。如图4所示。

流域的边界面由物面SB，尾涡面SW和外边界面S∞组成。在该流场中可用扰动速度势φ来表示螺旋桨的扰动。在oxyz坐标系中，Φ满足Laplace方程：

在螺旋桨表面满足流体法向速度为零的运动边界条件，即

式中：nQ是边界面上的单位法向量，在oxyz坐标系中表示。

图4 螺旋桨及周围流场示意图Fig.4 Propeller and fluid around it

1.3 时变尾涡数学模型

对于螺旋桨叶片这种升力体，从叶片随边泄露出的尾涡会影响叶片表面的环量，因此需要考虑尾涡的影响。对尾涡的建模主要考虑尾涡的强度和尾涡的形状。对于传统的根部固定不随轴系运动的螺旋桨，尾涡的形状通常假定为螺旋桨面，泄露的强度通常按Morino库塔条件或压力库塔条件处理[15]。而对于随轴系做复杂空间运动的螺旋桨，由于其运动轨迹复杂，使得不同时刻随边的位置也不一样，因而从随边泄露的第一个尾涡位置也不断变化。因此需要对尾涡进行合理建模，以便能考虑不同时刻螺旋桨位置变化对尾涡几何形状的影响。

在本文中，假设尾涡的泄露是一个按时间变化的过程，如图5所示。在初始时刻，假设螺旋桨静止，此时没有尾涡泄露。在Δt时刻，螺旋桨往前移动一个距离，此时泄出第一个尾涡，其强度用简单的Morino库塔条件，即：

式中：Δφ（rT）为叶片半径rT的随边处泄露的尾涡速度势，φ+（rT）为叶背（吸力面）随边处的速度势，φ-（rT）为叶面（压力面）随边处的速度势。

图5 尾涡泄露过程Fig.5 The process of vortex-shedding

在2Δt时刻，螺旋桨继续往前移动一个距离，此时泄出第二个尾涡，强度仍然按（6）式确定。而Δt时刻泄露的尾涡其强度保存不变，在原地运动并发生收缩、卷曲等变形。以此时间类推，从而尾涡的泄露是一个连续的过程。

泄露的尾涡是不受力的，由库塔-茹科夫斯基定理可知，泄露的尾涡的速度必定与当地的流场速度平行。也就是说泄露的尾涡片必定按当地的流线运动。在OXYZ坐标系中，尾涡的运动速度为Vw=Δφ。因此泄露的尾涡每个时间步运动的距离为：

2 求解过程

本文采用基于扰动速度势的面元法求解（4）式与（5）式组成的定解问题。根据Green定理，当场点P（x,y,z）在物面上时，扰动速度势可以表示为（结合边界条件）

式中，Δφ（Q1）为尾涡面的上速度势跳跃，按（6）式确定。

对于复杂形状结构，（8）式很难求解析解，因此采用数值面元法求解。假设螺旋桨共有Z个叶片，将一个叶片及相应轮毂划分成Np个四边形面元（沿弦向面元数为N，沿展向面元数为M，轮毂面元数为Nh，则Np=N·M+Nh）。泄露尾涡面元的展向数目为M，弦向数目由时间总步数Nt决定。比如在it个时间步，弦向泄露共it个尾涡面元，其中只有紧靠叶片随边的第一列面元强度未知，其它尾涡面元强度均已知。如图6所示。

图6 螺旋桨及泄露尾涡面元分布Fig.6 Panel arrangement of propeller and wake

因此在第it个时间步，（8）式可离散成：

由此可知，（9）式右边各项均为已知量，方程可以求解。进一步，（9）式写成矩阵形式：

式中：“±”表示随边处上下面元。

在oxyz坐标系中，不定常Bernoulli方程可表示为：

式中：P0为参考点处流体压力；ρ为流体密度；V＝Vin+▽φ，为总扰动速度。

求出叶片上速度势分布后，即可由（13）式求出叶片表面的压力分布，结合结构动力学，进而可以求解螺旋桨叶片变形、轴承力等等。整个求解流程如图7所示。

图7 求解流程图Fig.7 The process of solution

3 数值算例

本节针对实际使用中的螺旋桨给出两个算例，以验证算法的准确性，同时预报轴系振动对螺旋桨激励力的影响。从求解流程图7可以看出首先必须确定进流速度Vin，而由（3）式可知，需要知道轴系振动速度Vs与ωs，才能求出进流速度Vin。由于实际中轴的振动形式比较复杂，一方面受多种激励源影响（比如不平衡激励，弯曲激励，流体激励等等）；另一方面激励力与振动之间还存在耦合效应，比如流体激励力引起轴系振动，而轴系振动又对流体激励力有反馈作用。因此需要建立一个完整的流体-桨-轴系流固耦合动力学模型，才能较为准确地预报轴系振动，这已经超出本文的研究范围。为了研究方便，本文假设轴系的振动形式已知，为简单的简谐振动。本文以4381螺旋桨为研究对象，其为5叶片无侧斜螺旋桨，设其直径为5 m，其余详细几何参数在参考文献[16]中给出。

图8 脉动推力系数与扭矩系数比较Fig.8 The comparison of thrust and torque coefficients

第二个算例是研究螺旋桨在纵向与横向（垂向）联合振动下的脉动力变化。设轴系末端的振动速度为（单位为m/s）：

式中ω既是轴的振动频率，也是螺旋桨的旋转频率。因此可以认为轴系按转频振动。计算中设ω为3 Hz。假设轴系纵向振动幅值为1 mm，根据简谐振动理论，轴系纵向振动速度大约为0.018 m/s。

图9 脉动推力系数与扭矩系数Fig.9 The pulsant thrust and torque coefficients

图10 脉动推力与扭矩Fig.10 The pulsant thrust and torque

图11 脉动垂向力与横向力Fig.11 The pulsant vertical force and lateral force

螺旋桨转动一圈的过程中计算360步，一共计算6圈。图9为计算的推力系数与扭矩系数随旋转角度的变化过程。从图中可以看出，由于螺旋桨从静止启动，在启动瞬间，螺旋桨加速度无穷大，在桨叶表面瞬时产生一个启动涡，导致桨叶表面的脉动压力变化非常剧烈。在螺旋桨旋转大约两圈后，桨叶表面的压力脉动呈现周期性变化，表明此时计算已经趋于稳定。计算结果表明在（14）式的速度假设下，稳定后的脉动推力呈周期变化，变化频率为轴频，其推力脉动分量大约为静态分量的4.5/1 000，扭矩脉动分量大约为静态分量的4/1 000。对于大型船舶轴系来说，纵向振动1 mm的量级是非常普遍的，由此可见其产生的脉动力将对轴系振动产生重要影响，因而不能忽略。

图10为计算稳定后，叶片旋转一圈的脉动推力与扭矩变化。图11为脉动的垂向力与侧向力变化。此处，脉动力均去掉了静态分量，只考虑其脉动部分。从图中可以看出，单个叶片上旋转一圈的载荷不是按简谐变化的，表明此时有轴频以外的其它频率成分存在，进一步研究后表明是轴频的倍频成分。此频率的来源是（13）式中计算脉动压力时，对速度进行平方而引入的。但是所有叶片的合力又是一个简谐激励力，表明不同叶片的倍频分量相位差满足严格的滞后关系，导致它们可以相互抵消。又由于所有叶片脉动力的轴频分量的相位均一致，使得总的合力有所增强。

4 结 语

本文构建了螺旋桨在流场中随轴系做复杂空间运动时的水动力预报数学模型，并利用不定常面元法求解了这一问题。通过与已有文献的研究结果做比较，证明了本文算法的有效性。利用本文的方法可以方便地研究螺旋桨随轴系振动而做复杂空间运动时的水动性能，分析轴系振动对螺旋桨水动力的影响，从而解决了上文提到的流体-螺旋桨-轴系流固耦合问题中的第三个关键子问题，为研究流体-螺旋桨-轴系双向流固耦合动力学特性做基础。同时该方法也可以用来预报船舶转弯、船舶升沉及纵摇振荡、螺旋桨启停及加速等复杂工况下的螺旋桨水动性能。最后预报了轴系纵向振动幅值1 mm，振动频率为3 Hz并伴有小量回旋振动下4381螺旋桨的脉动力，结果表明推力脉动分量大约为其静态分量的4.5/1 000，扭矩脉动分量大约为其静态分量的4/1 000。因此轴系振动对螺旋桨脉动力的影响不能忽略。