1. 引言

直翼推进器是一种回转轴垂直设置的新型推进装置 [1] ，其良好的操纵性及较优的推进效率使之广泛适用于现代化拖船 [2]、海洋深潜器 [3]、港作清洁船 [4]、车客渡船、游览船等民用船舶。针对直翼桨型船舶开展振动特性分析具有一定的工程实用价值，然而目前还没有针对直翼桨的激振力经验公式。

目前针对常规推进器激励下的尾部振动研究较为丰富 [5] ，大多数学者利用经验公式来确定螺旋桨诱导的船底表面力。刘西安 [6] 基于国军标得到了螺旋桨产生的表面力，将其作为输入载荷，计算并评估了该船尾部结构的振动响应。周清华等 [7] 通过经验公式得到了螺旋桨脉动压力载荷，采用直接瞬态响应分析方法结合虚拟质量法研究了某滑行艇尾部结构的动力响应。朱胜昌等 [8] 根据DNV公式确定了五、六缸机的螺旋桨激振力，使用模态叠加法对上层建筑开展了强迫振动计算。由于昂贵的时间和经济成本，少数学者根据实测数据来预报船尾结构的振动响应，许树浩等 [9] 利用实船测试得到测点的脉动压力及各点之间的相位差，并以此作为输入进行了振动响应计算。朱理等 [10] 使用内、外桨区试验得到的脉动压力结果对某舰进行了振动特性分析。迄今为止，德国、美国、日本等国的科研人员已对直翼推进器展开了大量的水动力性能研究，并取得了许多重要成果 [11] [12] [13] [14] [15]。 虽然上述对常规推进器激励下的尾部振动研究及直翼桨水动力性能研究的文献较多，但对直翼推进器激励下的船体尾部振动研究却很少。

本文以直翼桨激振力作为输入载荷，利用结构有限元法结合虚拟质量法，计算船舶尾部结构在额定工况下的振动响应，分别使用频域响应分析方法和时域分析方法计算尾部振动，验证频域响应分析方法的准确性。为了减小尾部振动响应，由初始船型出发，从改变板厚、加强基座结构及修改中纵舱壁等可行性方法进行优化，对比各方案在额定工况下的速度响应。

2. 尾部结构振动响应的模态叠加法

通常由推进器作用在船舶尾部结构的载荷为随时间变化的载荷，求解振动响应可以在时域或频域内进行。模态叠加法是一种求解动力响应的数值方法，其利用结构的模态振型来对耦合的运动方程进行缩减和解耦，然后再由单个模态响应的叠加得到问题的最终解答。对于多工况的复杂模型而言，模态法在缩短计算周期的同时又能够达到较高的计算精度。由于船舶尾部结构较为复杂，对其振动响应求解采用模态叠加法。

2.1. 频域响应分析理论

简谐激励下有阻尼多自由度系统在模态坐标下的运动方程为 [16]：

$-{\omega }^2\left[M\right]\left[\varphi \right]\left\{\xi \left(\omega \right)\right\}+i\omega \left[B\right]\left[\varphi \right]\left\{\xi \left(\omega \right)\right\}+\left[K\right]\left[\varphi \right]\left\{\xi \left(\omega \right)\right\}=\left\{p\left(\omega \right)\right\}$(1)

式中： $\omega$为激励力频率， $\left[\varphi \right]$为坐标变换矩阵， $\left\{\xi \left(\omega \right)\right\}$为模态坐标向量， $\left[M\right]$为质量矩阵， $\left[B\right]$为阻尼矩阵， $\left[K\right]$为刚度矩阵， $\left\{p\left(\omega \right)\right\}$为激励力幅值向量。

前乘，将方程解耦后得到：

$-{\omega }^2{\left[\varphi \right]}^{\text{T}}\left[M\right]\left[\varphi \right]\left\{\xi \left(\omega \right)\right\}+i\omega {\left[\varphi \right]}^{\text{T}}\left[B\right]\left[\varphi \right]\left\{\xi \left(\omega \right)\right\}+{\left[\varphi \right]}^{\text{T}}\left[K\right]\left[\varphi \right]\left\{\xi \left(\omega \right)\right\}={\left[\varphi \right]}^{\text{T}}\left\{p\left(\omega \right)\right\}$(2)

式中： ${\left[\varphi \right]}^{\text{T}}\left[M\right]\left[\varphi \right]$为模态质量矩阵， ${\left[\varphi \right]}^{\text{T}}\left[B\right]\left[\varphi \right]$为模态阻尼矩阵， ${\left[\varphi \right]}^{\text{T}}\left[K\right]\left[\varphi \right]$为模态刚度矩阵， ${\left[\varphi \right]}^{\text{T}}\left\{p\left(\omega \right)\right\}$为模态力向量。

当使用模态阻尼时，运动方程保持为非耦合形式，每一模态有如下形式方程：

$-{\omega }^2{m}_i{\xi }_i\left(\omega \right)+i\omega b_i{\xi }_i\left(\omega \right)+k_i{\xi }_i\left(\omega \right)=p_i\left(\omega \right)$(3)

式中：下标i指第i阶模态。

其中每一个模态的响应由下式计算：

${\xi }_i\left(\omega \right)=\frac{p_i\left(\omega \right)}{-m_i{\omega }^2+i{b}_i\omega +k_i}$。(4)

2.2. 时域响应分析理论

通过使用中心差分法计算离散点处的速度 $\left\{{\stackrel{\dot{}}{u}}_n\right\}$及加速度 $\left\{{\ddot{u}}_n\right\}$，外力取三个邻近时间点的平均值，将变量从物理坐标转换到模态坐标后，得到：

$\left[A_1\right]\left\{{\xi }_{n+1}\right\}=\left[A_2\right]+\left[A_3\right]\left\{{\xi }_n\right\}+\left[A_4\right]\left\{{\xi }_{n-1}\right\}$(5)

$\{\begin{array}{l}\left[A_1\right]=\left[{\varphi }^{\text{T}}\right]\left[\frac{M}{\Delta t^2}+\frac{B}{2t}+\frac{K}{3}\right]\left[\varphi \right]\\ \left[A_2\right]=\frac{1}{3}\left[{\varphi }^{\text{T}}\right]\left\{P_{n+1}+P_n+P_{n-1}\right\}\\ \left[A_3\right]=\left[{\varphi }^{\text{T}}\right]\left[\frac{2M}{\Delta t^2}-\frac{K}{3}\right]\left[\varphi \right]\\ \left[A_4\right]=\left[{\varphi }^{\text{T}}\right]\left[-\frac{M}{\Delta t^2}+\frac{B}{2t}-\frac{K}{3}\right]\left[\varphi \right]\end{array}$

式中： $\left[A_1\right]$为动力矩阵， $\left[A_2\right]$为外力矩阵， $\Delta t$为积分时间步长。

当考虑模态阻尼时，运动方程保持为非耦合形式，每一模态有如下形式方程：

$m_i{\xi }_i\left(t\right)+b_i{\xi }_i\left(t\right)+k_i{\xi }_i\left(t\right)=p_i\left(t\right)$(6)

其中每一个模态的响应由通解及特解组成：

${\xi }_i\left(t\right)={\text{e}}^{-b_{i}t/2m_i}\left({\xi }_0\cos {\omega }_{d}t+\frac{{\xi }_0+b_i/\left(2m_i\right){\xi }_0}{{\omega }_d}\sin {\omega }_{d}t\right)+{\text{e}}^{-b_{i}t/2m_i}\frac{1}{m_i{\omega }_d}{\displaystyle {\int }_0^t{\text{e}}^{b_i\tau /2m_i}}p\left(\tau \right)\sin {\omega }_d\left(t-\tau \right)\text{d}\tau$。(7)

3. 有限元模型

在采用有限元法对船舶尾部振动问题开展研究时，通常对船体结构使用全船模型、尾部舱段模型及中、首船体梁模型与尾部模型相结合的整船混合有限元模型。本文研究对象为某双直翼桨船舶局部区域，由于尾部局部模态的计算受整体结构有限元模型差异影响较小 [17] [18] ，因此可直接使用合适的尾部舱段模型进行计算。该模型相较于全船三维空间有限元模型或混合模型，建模效率更高、计算耗时更短，同时也能够保证计算结果的准确性。本文参照规范 [19] 中的相关规定，利用MSC. Patran有限元软件对尾部舱段模型进行建模，三维有限元模型见图1。

尾部舱段模型共有17,341个单元，结构总长L = 30.6 m，型宽B = 6.8 m，型深D = 2.1 m。纵向范围为尾封板向船艏延伸至1/4倍船长，包含3个舱室的长度，横向及垂向范围采用全宽、全高模型。模型X轴从船艏指向船尾，Y轴指向右舷，Z轴指向甲板，坐标原点在模型最前端肋位处基线位置。将模型最前端剖面上节点的位移边界条件设为简支。

单元类型主要为四节点四边形板壳单元、梁单元以及小部分的三角形单元。外板、舱壁、甲板、平台、基座和强框架等结构使用板壳单元模拟，其余骨材均用梁单元模拟。通过对单元材质密度的修改保证有限元模型与实船质量和分布一致。船上主要设备及直翼螺旋桨的重量使用质量单元模拟，将质量单元与甲板上设备安装位置对应的节点采用MPC多点约束的方法进行连接，设置MPC类型为RBE2刚性单元。附连水质量采用虚拟质量法，通过定义湿表面单元和吃水高度来实现耦合振动计算 [20]。

4. 直翼推进器激振力

直翼推进器又称摆线推进器，叶片与传统推进器的螺旋叶片不同，呈剖面为翼型的直叶片，一般由4至6片桨叶组成，其伸出船体表面并与船底相垂直，即轴的方向为垂直方向，故又称竖轴推进器。船舶行进时，其桨叶绕叶轮垂直轴线匀速旋转，同时每片桨叶分别绕着自身旋转轴做角摆动，而叶轮又随船体做直线运动，三者的合运动为摆线运动 [21]。 在实际应用中一般安装两个相向旋转的直翼推进器，以达到转矩互相抵消的目的，同时利用两个推进器的推力方向控制，操控船舶运动，使之转向、横移和原地回转。其实物及实际安装情况如图2所示。

推进器的激振力由轴承力和表面力组成，沿船体表面对脉动压力积分得到的推进器对船体总激振力称为表面力，通过轴系和轴承传递给船体，作用于桨叶上的变动流体力所引起的激振力成为轴承力，其中因空化现象产生的脉动压力是引起船舶尾部结构振动的主要振源 [23]。 目前确定激振力的主要途径有理论计算、经验公式、模型试验和数值模拟。近年来，虽然学者们针对螺旋桨空泡特性方面开展了大量的研究，但由于船、桨、舵之间复杂的耦合作用，理论计算方法尚存在一定缺陷。经验公式由大量试验数据归纳总结得到，是计算激振力较为常用的一种手段，但其缺点在于计算精度不够。模型试验和CFD数值计算是学者们广泛认可的两种方法。常规试验时，螺旋桨诱导脉动压力的测量在空泡水筒中进行，以船模后半体和金属网格来模拟伴流场，一般在螺旋桨上方平板上布置4至5个压力传感器，或者在螺旋桨上方和前方的船模后半体上布置3至4个压力传感器 [24]。 相较于数值模拟方法，模型试验的研究成本较高、研究周期较长，但其能够模拟实船伴流场中螺旋桨的工作状态，得到高精度的数据结果。

4.1. 直翼推进器轴承力

直翼推进系统包括直翼推进器、传动轴、齿轮箱和主机等。由于直翼推进器在非均匀的进流中旋转，导致桨叶上推力和旋转阻力变化，因此产生了脉动分力和转矩。由于直翼桨与船体之间使用螺栓连接，在本模型中利用MPC多点约束来模拟脉动推力及转矩的传递，轴承力作用点在直翼桨的旋转轴线上，沿船底向下延伸1/2直翼桨叶片长度，桨的有限元模型见图3。

以一阶叶频为例，本文根据模型试验提供的数据绘制了4叶型直翼桨轴承力幅值与转速之间的关系图，试验在如图4所示。由图可见，三个脉动分力及脉动转矩幅值随转速的增大而增大。

4.2. 直翼推进器表面力

直翼推进器旋转时产生压力场，压力场通过其周围的流体介质传递作用到船体表面。经过对实船或船模的测试资料进行统计分析，专家们总结出一些可用于初步预估传统螺旋桨脉动压力的近似算法，例如挪威船级社的霍尔顿方法、日本学者谷口中——高桥肇方法、藤野方法、O.M.帕利等编写的《船舶结构力学手册》中的方法等 [25] ，然而这些经验公式并不适用于直翼推进器。以右舷为例，本文根据设备厂商提供的试验数据，将直翼桨表面力作用范围划分为15个区域，分布如图5所示，T1至T4为直翼桨前端区域、T5至T7为直翼桨外侧区域、T8至T10为直翼桨内侧区域、T11至T13为船中区域、T14至T15为直翼桨后端区域，试验在空泡水筒中进行。

为了模拟实际压力分布情况并综合考虑船底板架网格形式，将表面力加载到其对船底相应的作用区域，如图6所示。

5. 尾部振动响应计算分析

对尾部三维有限元模型分别使用模态频域响应分析和时域响应方法计算其在额定工况下直翼桨激振力作用下的振动响应，并考虑不同载荷间相位角的影响。参照规范 [26] 中船舶振动测量方法，选取甲板中纵剖面尾端处为典型节点，提取该节点在两种方法下的计算结果。当进行时域响应分析时，应当根据激励载荷的大小选择合适的时间步，步长和采样点数，保证每个周期内至少有8个时间步 [27] ，随着步长的减小，求解精度逐渐提高，但同时也会导致计算时间增长，而时间步及步长过短会导致响应未达到稳态。根据试算结果，选取时间步为3000步，步长为0.005 s。将时域响应分析的结果进行快速傅里叶变换，并与频域响应分析得到的数据进行对比，得到两种计算方法的相对误差，计算结果见表1。

由表1的计算结果可以看出，使用频域响应分析方法及时域响应分析方法计算得到的整体相对误差较小，在两种叶频工况下的误差不超过±3%，结果较为可靠。

6. 尾部结构多方案振动响应对比

为了研究尾部结构对该船型振动响应的影响，本文基于原始方案，通过改变板厚、加强基座结构及修改中纵舱壁等可行性方法进行优化，使用模态频域响应分析法计算各方案在额定工况下的速度响应。方案一为原始方案，其余方案均基于方案一进行优化；方案二将直翼桨作用区域至尾部结构末端之间水线以下船体外板板厚增加4 mm，尾部结构总重量增加了5.47%；方案三加强了直翼推进器基座部分结构，将每个基座四周设置了四根高腹板加强筋，尾部结构总重量增加了0.627%；方案四在去除直翼推进器基座区域中纵舱壁结构的同时设置支柱，尾部结构总重量减轻了0.336%，计算结果见表2、表3。

由表2、表3的计算结果可以看出：

1) 参考规范 [26] ，方案一至方案四速度响应幅值最大点的垂向、横向或纵向的振动参数均低于振动衡准值7 mm/s，满足规范要求，故四个方案均为有效方案，其中方案四的尾部结构总重量最小，原始方案、方案三次之，方案二最大。各方案的垂向振动响应最为显著，总体而言，方案四的速度响应幅值最大点垂向速度响应分量相对较大，方案二的速度响应幅值最大点垂向速度响应分量相对较小。

2) 引起该船型尾部振动的主要激励力是直翼桨叶频分量，因此以各方案在叶频激励下的振动响应为评判减振效果的主要依据。由方案二可以看出，增加船体底板板厚对于减小尾部振动响应具有一定的效果，在激励频率为叶频及倍叶频时速度响应幅值最大点的垂向速度响应分量分别降低了4.79%、12.61%；由方案三可以看出，增强直翼推进器基座部分结构对于减小尾部振动响应同样具有一定的效果，在激励频率为叶频时速度响应幅值最大点的垂向速度响应分量降低了1.68%。对比方案一和方案四可知，在激励频率为叶频时速度响应幅值最大点的垂向速度响应分量相差了0.062 mm/s，纵舱壁的设置对尾部结构固有频率具有较大的影响，对纵舱壁进行合理设计可大幅度降低振动响应。因此，在综合考虑结构总重量及振动响应情况的前提下，建议采用方案四。

7. 结论

本文针对直翼桨激励下的船舶尾部振动问题，采用尾部三维有限元模型，使用基于模态叠加法的频域响应分析和时域响应分析，计算额定工况下的振动响应，并考虑载荷间相位角的影响。根据计算结果可知两种分析方法之间的相对误差最大约±3%，计算结果较为一致。基于原始模型进行多方案设计，分析不同尾部结构方案下的振动影响，由计算结果可知引起船体尾部振动的主要激励力是叶频分量，合理设置纵舱壁可大幅度降低振动响应，适当增加板厚也具有一定的尾部减振作用。

基金项目

国家自然科学基金资助项目(51839005，51879113)，中央高校基金科研业务费资助项目(2019kfyXMBZ048)。

NOTES

^*通讯作者。

参考文献