基于Newmark-β法的鱼雷推进轴系动力响应分析

赵 琪，尹韶平，王 中，郭 君

(1．中国船舶重工集团公司第705所,西安 710075；2.水下信息与控制重点实验室，西安 710075)

基于Newmark-β法的鱼雷推进轴系动力响应分析

赵 琪1，2，尹韶平1，2，王 中1，2，郭 君1

(1．中国船舶重工集团公司第705所,西安 710075；2.水下信息与控制重点实验室，西安 710075)

针对轴系在不同激励条件下的动力响应问题，基于Newmark-β法，运用MATLAB软件编程，实现对鱼雷推进轴系的动力响应分析，通过与ANSYS软件计算结果的对比，验证所建模型和MATLAB计算程序的合理性；同时得出推进轴系的动力响应与激励条件的线性变化关系。

鱼雷推进轴系;动力响应;Newmark-β法

当轴系受到动载荷的作用时会产生振动，引起轴系局部疲劳，甚至使轴系破坏，因此必须进行动力分析。应用有限元理论对轴系进行动力分析，一般分为振动特性分析和动力响应分析两种。轴系的振动特性反映轴系的固有特性，是研究动力问题的基础，一般包括频率和振型的计算。动力响应分析是指轴系在外力作用下的强迫振动，主要求解轴系的位移、速度、加速度等物理量随时间的变化情况[1]。

文献[2]应用MATLAB软件建立了鱼雷推进轴系的理论模型，对模型进行了振动特性分析，并通过与ANSYS软件计算结果的对比，验证了所建模型的合理性。在此基础上文中采用Newmark-β法对鱼雷推进轴系进行动力响应分析，一方面通过与ANSYS计算结果进行对比，对模型进一步验证；另一方通过分析鱼雷推进轴系的动力响应与输入参数的变化关系，对轴系设计提出要求。

1 鱼雷推进轴系理论模型

鱼雷推进轴系结构见图1。

图1 鱼雷推进轴系结构

动力装置与花键轴之间、花键轴与尾轴之间通过花键相连，推进器转子安装在尾轴上。花键轴及尾轴主要用于将动力装置的输出功率传递给推进器，而推进器工作时产生的推力通过推力轴承传递至雷尾壳体，推动鱼雷航行。在花键轴与尾轴相连的一端设有弹性花键套，将刚性连接转换为柔性连接。尾轴的支撑是靠隔板的滚珠轴承和尾盖处的滑动轴承，在滚珠轴承及滑动轴承外侧都安装了金属橡胶隔振器。

按照转子动力学分析的有限元方法，建立鱼雷推进轴系的理论模型，见图2。

图2 鱼雷推进轴系理论模型

将鱼雷推进轴(花键轴和尾轴)的质量、转动惯量离散化到结点上；联轴器利用等效轴段法进行处理；用刚度阻尼系数矩阵模拟轴承与隔振系统；轴系支撑结构通过等效刚度与等效质量单元代替。从左至右共划分22个结点，取轴段变截面处及支撑点位置。结点1为花键轴前端；结点8、9、10为弹性花键套简化位置；结点11、19为尾轴支撑位置。

基于此模型，建立系统运动微分方程[3]。

(1)

式中：M、G、K分别为系统整体质量矩阵、回转矩阵和刚度矩阵，由组成系统的各单元矩阵按照对应结点叠加而成[4]；Q为干扰力即激振力向量；U为结点位移列阵。

其中：回转矩阵G为转速Ω与整体极转动惯量矩阵J的乘积，即G=ΩJ。

整体质量矩阵、回转矩阵和刚度矩阵，在不同的振动形式下会有所不同[5-6]。文献[2]按照轴系在纵向、扭转、弯曲耦合振动状态(即单元结点为六自由度)下，确定了各自的矩阵形式。文中将利用所建模型，基于Newmark-β法对系统进行动力响应分析，即研究鱼雷推进轴系轴系在不同激振情况下的振动特性。

2 Newmark-β法

2.1 Newmark-β算法机理

求解结构动力响应问题，常用的方法有振型叠加法和逐步积分法。振型叠加法是利用结构自由振动振型互相正交的特性，将结构动力学方程化成各广义坐标的非耦合方程，然后进行求解[7]。振型叠加法适用于像地震等只激发起较少振型，所需计算的响应历程较长这类问题。逐步积分法包括传统的中心差分方法、Wilson-θ法、Newmark-β法等。逐步积分法对于线性系统和非线性系统都适用，而且在MATLAB程序上更容易实现，所以文中采用Newmark-β逐步积分法计算结构的动力响应。

(2)

Newmark-β法假设t+Δt时刻的速度和位移存在如下关系式：

(3)

(4)

式中:0≤γ≤10≤2β≤1。

参数γ和β的选择对算法影响很大。算法稳定性分析指出，当γ≥0.5，β≥(1/2+γ)2/4时,Newmark-β法无条件稳定，这时可以只根据精度的要求选择时间步长Δt。

利用Newmark-β法进行计算，即可得到单元结点的动力响应，包括位移、速度和加速度值。

2.2 时间步长的选择

时间步长Δt的选择体现了算法的精度要求[8]。Δt应能够根据所输入的状态自动进行调整。假设激励频率为f(Hz)，则时间步长可设置为Δt=1/nf(n=1,2,…) ，总时间可设置为tend=m/f(m=1,2,…)。其中：n值的意义是在每个周期内选取进行模拟的点数，即n值越大，计算精度越高，但同时计算量也会增大，延长计算时间，所以要根据实际需要进行选取。

这种简便的动态时间步长，可根据所输入的状态自动进行调整，在编程和计算中更加简便和灵活。

3 轴系的动力响应分析

结合鱼雷推进轴系的实际运转情况，具体分析其在3种不同激励下的位移响应。

1)轴系动平衡状态。模拟轴系由于加工等原因造成的质量分布不均而产生的激励力。

2)花键联轴器不对中状态。将花键轴简化为动力装置与尾轴之间的联轴器。由于制造精度或安装误差等原因，使得花键轴与尾轴轴线不重合，产生不对中啮合力，对轴系振动特性产生影响。

3)强迫位移响应。模拟动力装置的跳动量对推进轴系的激励。

3.1 动平衡响应分析

鱼雷花键轴和尾轴由于原材料质量分布不均或加工精度等原因，会产生一定的偏心，这对推进轴系的振动特性会产生较大影响。这种动平衡状态在理论模型中可用结点处圆盘的偏心质量来表示，见图3。圆盘的重心c与转轴的中心o′不重合。

图3 圆盘的偏心

当圆盘以角速度ω转动时，重心c的加速度在坐标轴上的投影为

(5)

式中：e——圆盘的偏心距，e=o′c。

在转轴的弹性力F作用下，由质心运动定理，有

(6)

带入式(5)得轴心o′的运动微分方程

(7)

这为由轴系偏心质量引起的强迫振动的微分方程。

假设偏心质量为m0，相位角φ0，则由于偏心质量即不平衡质量引起的激振力表达式为

(8)

将激励力带入系统运动方程，即可求解结点处位移响应。具体计算结果如下。

选取靠近花键轴中部(结点5)、靠近尾轴中部(结点14)在变转速(0～500Hz)下的位移分析曲线。图4a)为MATLAB计算结果，图4b)为ANSYS计算结果。

图4 动平衡响应

由图4可见，MATLAB曲线峰值对应的频率点较ANSYS偏大，但总体趋势较为接近。

动平衡响应可控因素有不平衡质量和偏心距，下面分析位移响应与这两种因素的变化关系。见图5。

图5 位移响应

由图5a)可见，结点处位移响应随不平衡质量的增大而线性增大。

由图b)可见，结点处位移响应随偏心距的增大而线性增大。

3.2 不对中响应分析

将花键轴简化为动力装置与尾轴之间的联轴器，联轴器的平行不对中见图6。

图6 花键联轴器的平行不对中

图7 不对中状态的花键啮合情况

由于花键的不对中，会使得各花键在转动过程中产生变形，影响啮合刚度。在联轴器传递转矩的过程中，由于扭转以及两个半联轴器的动态振动位移，都会产生啮合力，这都会对鱼雷推进轴系的振动特性产生影响。

此时，联轴器不对中啮合力fx、fy为

(9)

Fx、Fy及φ的具体推导过程参考文献[9]。

对于花键联轴器，当联轴器存在不对中时，轴套会产生偏心，在系统运行过程中，产生不平衡激励，且该激励的频率为转速频率的2倍，轴套偏心激励表达式如下。

(10)

式中：2m——花键轴的质量。

因此，花键联轴器不对中状态下的广义力应包括不对中啮合力和轴套的不平衡激励力。

花键联轴器不对中啮合力与静态不对中量和动态不对中量有关。假设无动态位移，则不对中啮合力随静态位移的变化关系见图8。

图8 啮合力随静态位移的变化

由图8可见，啮合力随静态位移的增大呈线性增加。

联轴器动态不对中量可控制程度较小，所以可分析位移响应与静态不对中量之间的变化关系，见图9。

图9 位移响应与静态不对中量的变化关系

图9中静态不对中量e0在0～2.5 mm之间的变化。由图9可见，不对中位移响应随静态位移的增大而线性增大。

3.3 强迫位移响应分析

与鱼雷尾舱段相连的动力舱段中，对动力装置采取了整机隔振的措施，使得动力装置产生了较大的跳动量。这种跳动给推进轴系的振动特性带来较大的影响，它的作用可以利用强迫位移模型进行分析。

假设强迫位移幅值为A，则在x、y方向位移可写为

(11)

由位移表达式可得x、y方向上的加速度，所以由强迫位移引起的激励力为

(12)

将激励力带入系统运动方程，即可求解结点处位移响应。具体计算结果如下。

选取靠近花键轴中部(结点5)、靠近尾轴中部(结点14)，在变转速(0～500 Hz)下的位移分析曲线。MATLAB计算结果见图10，ANSYS计算结果见图11。

图10 强迫位移响应MATLAB计算结果

图11 强迫位移响应ANSYS计算结果

由图10和图11可见，两者曲线较为接近，且峰值对应的频率点也基本相同。

轴系位移响应随动力装置跳动量的变化情况见图12。

图12 位移响应与强迫位移的变化关系

由图12可见，强迫位移幅值A在0～5 mm之间变化，结点处的位移响应与强迫位移的变化是线性变化的，并且，越靠近激励点，其位移响应越大，且变化较为明显。为了控制轴系的强迫位移响应，可对动力装置的跳动量提出具体的要求。

4 结论

文中利用所建鱼雷推进轴系理论模型，基于Newmark-β法对其进行动力响应分析，通过与ANSYS软件计算结果进行对比，进一步验证了模型的正确性。基于此模型，可以对轴系其它振动特性进行分析，也可对模型的参数进行优化处理，对轴系初始设计提出指导意见。

此外，通过分析轴系激励可控参数可以看出，轴系结点的位移响应与不平衡质量、偏心距、联轴器静态不对中量、强迫位移量等可控参数之间，呈现出线性变化关系。根据这种线性变化关系，一方面可针对轴系结构设计合理的减振降噪装置；另一方面在轴系的初始设计阶段，可对联轴器静态不对中量、动力装置的跳动量等提出合理的设计要求。

[1] 谢龙汉,刘新让,刘文超.ANSYS结构及动力学分析[M].北京:电子工业出版社,2012.

[2] 赵 琪,尹韶平,王 中,等.基于MATLAB的鱼雷推进轴系弯曲振动涡动频率计算[J].鱼雷技术,2014,22(6):221-301-308.

[3] 钟一谔.转子动力学[M].北京:清华大学出版社,1984.

[4] 王勖成.有限单元法[M].北京;清华大学出版社,2003.

[5] 陈锡恩,高 景.船舶轴系回旋振动计算及其参数研究[J].船海工程,2011(5):8-11.

[6] 岳 聪,任兴民,邓旺群.柔性转子加速过临界瞬态响应特征分析[J].机械科学与技术,2013,32(3):395-398.

[7] 李鸿晶,王 通,廖 旭.关于Newmark-β法机理的一种解释[J].地震工程与工程震动,2011,31(2):55-62.

[8] 谢国伟,黄亚宇.一种自适应步长的Newmark-β预测校正积分法[J].新技术新工艺,2013(5):6-8.

[9] 赵 广.转子-联轴器-轴承-隔振器系统耦合动力学特性研究[D].哈尔滨:哈尔滨工业大学,2009.

Analysis of Dynamic Response for the Torpedo Propulsion-shaftSystem Based on Newmark-β Method

ZHAO Qi1,2, YIN Shao-ping1,2, WANG Zhong1,2, GUO Jun1

(1 The 705 Research Institute, China Shipbuilding Industry Corporation, Xi′an 710075, China;2 Science and Technology on Underwater Information and Control Laboratory, Xi′an,710075, China)

In order to study dynamic response of the shafting based on multiple-excitation, the MATLAB program is used to achieve the response analysis of torpedo propulsion shafting on the basis of Newmark-β method. By comparing the results with that of ANSYS, the reasonableness of the model and the MATLAB program is verified. The linear relationship between analysis of dynamic response and input conditions is also given.

Newmark-β; torpedo propulsion-shafting; dynamic response

10.3963/j.issn.1671-7953.2015.01.030

2014-09-11

国家部委重大专项

赵 琪(1989-)，男，硕士生

U661.31

A

1671-7953(2015)01-0117-05

修回日期：2014-09-25

研究方向:鱼雷总体技术

E-mail:qz_zhao@163.com