弧齿锥齿轮系统的非线性动力学特征分析*

曲凯明,郭小虎,谢振东,丁顺良,高建设

(1.郑州大学机械与动力工程学院,郑州 450001;2.格力电器(郑州)有限公司,郑州 450001;3.威海重鼎托邦智能装备有限公司,威海 264200)

0 引言

齿轮系统是机械装备中应用最广泛的动力和运动传递装置,其动力学行为和工作性能对整个机械装备有至关重要的影响[1]。齿轮系统动力学是研究齿轮系统在传递动力和运动过程中动力学行为的一门学科,在齿轮研究领域一直受到科研人员广泛关注[2-3]。弧齿锥齿轮是重要的基础传动部件,与直齿轮相比于,其具有承载能力强、传动比高、工作可靠和传动平稳等优点,目前广泛应用于汽车、航空航天、精密机械等领域[4-5],因此研究弧齿锥齿轮的非线性动力学特征十分重要。

目前对直齿轮、行星轮的非线性动力学特征已经有了较成熟研究[6-8]。而弧齿锥齿轮由于模型较为复杂,目前对其非线性动力学特征的研究相对还不够系统。王三民等[9]推导了考虑啮合刚度和齿侧间隙的弧齿锥齿轮传动系统八自由度非线性振动模型。杨先勇等[10]利用庞加莱映射、相空间重构和最大Lyapunov指数等方法,对含间隙、时变啮合刚度和传动误差的弧齿锥齿轮系统的非线性动力学特征进行研究。冯刚等[11]利用ANSYS分析了锥齿轮有裂纹和无裂纹两种状态下的单齿啮合刚度,并通过FFT频谱图和庞加莱映射研究了裂纹对锥齿轮系统振动的影响规律。XU等[12]建立了考虑静态传动误差的锥齿轮系统八自由度数学模型,研究了啮合频率、负载系数和齿侧间隙对锥齿轮系统振动特性的影响。CHEN等[13]建立了具有挤压油膜阻尼器的锥齿轮系统七自由度动力学模型,并通过相空间重构和Lyapunov指数等方法证明了挤压油膜阻尼器的应用可以有效抑制跳跃现象。WANG等[14]建立了瞬时混合润滑条件下考虑薄膜润滑和粗糙接触的弧齿锥齿轮数学模型,研究了粗糙度、润滑、施加力矩等对锥齿轮系统非线性动力学特性的影响。

弧齿锥齿轮系统是高维复杂非线性系统,目前在其非线性动力学特征分析方面主要用到相空间重构、最大Lyapunov指数、庞加莱映射等方法,利用递归图和递归量化分析方法的较少[15-16],同时缺乏啮合频率、刚度系数和阻尼系数等重要影响因素的综合分析。而递归图和递归量化分析方法,能够更加直观、清晰地反应弧齿锥齿轮系统的运行状况,从而定性、定量地分析系统的非线性动力学特征。

本文首先通过弧齿锥齿轮系统的数学模型,得到系统在不同啮合频率、刚度系数和阻尼系数下的时间序列,然后利用递归图和递归量化分析方法对弧齿锥齿轮系统的在不同参数下的非线性动力学特征进行研究,并利用相空间重构方法进行验证。

1 弧齿锥齿轮系统的动力学模型

目前对弧齿锥齿轮的建模已经有了较多的研究。本文使用的弧齿锥齿轮系统动力学模型如图1所示[11]。

图1 弧齿锥齿轮系统示意图

在该模型中,两齿轮均视为刚体。齿轮1为主动轮,受到的驱动力矩为T1,T1由不变部分T1m和变化部分T1v两部分组成;齿轮2为从动轮,受到的力矩为T2。设两齿轮绕其中心点在3个坐标方向上的移动和转动分别为{X1,Y1,Z1,θ1,X2,Y2,Z2,θ2}T,两弧齿锥齿轮齿面啮合点因振动和误差而产生的沿啮合点法线方向的相对位移为[9]:

λn=(X1-X2)c1-(Y1-Y2)c2-
(Z1-Z2+r1θ1-rrθ2)c3-en(t)

(1)

式中:r1和r2分别为两齿轮啮合点的半径,en(t)为齿轮副的法向静态传动误差,en(t)的计算公式为:

(2)

式中:Ωh为啮合频率,Ael和Φel分别为误差的l阶谐波幅值和初始相位。设齿轮副在啮合时的法向动载荷沿其坐标方向上的分力分别为:

(3)

式中:ch为啮合阻尼,kh(t)为时变啮合刚度,f(λn)为间隙函数。kh(t)和f(λn)的计算公式分别为:

(4)

(5)

式中:km为时变啮合刚度的平均值,Akl和Φkl分别为误差的l阶谐波幅值和初始相位,bm为法向平均齿侧间隙的一半。其计算公式分别为:

式中:ΩF为外载激励角频率,AFl为外载荷的l阶谐波幅值,ΦFl为初始相位。图1得到的无量纲方程为:

(6)

式中:各无量纲的参数的计算公式如下:

c1=cosδ1sinαn
c2=cosδ1cosαnsinβm
c3=cosαncosβm
c4=sinαncosβm+cosαnsinβmsinδ1
c5=sinαnsinδ1-cosαnsinβmcosδ1
c6=cosαnsinβm

i=x,y,z,j=1,2。时变刚度和误差激励的高阶分量对系统影响较小,本文取基频分量进行计算,即:

式中,ε和γ分别为刚度系数和静态传动误差系数。

2 非线性分析方法

2.1 相空间重构

相空间重构是一种用于分析复杂系统动力学特征的数学方法,其认为系统中每个分量的演化过程都包含着系统的全部信息。相空间重构方法能够在保持原始序列动态行为的前提下,挖掘时间序列所隐含的信息,适于处理较长的非平稳时间序列。根据Takens嵌入定理[17],一维离散时间序列{x(t),t=1,2,…,n},做变换后得到:

Xt=[x(t),x(t+τ),x(t+2τ),…,x(t+(m-1)τ)]
t=1,2,…,N

(7)

式中:τ为延迟时间,m为嵌入维数。根据相空间重构原理,可以将时间序列转换成:

(8)

从上式中可以看出,对于N个数据点的时间序列,若延迟时间为τ,嵌入维数为m,经重构后可得到N-(m-1)τ个向量。相空间重构需要确定两个重要参数:延迟时间τ和嵌入维数m。

2.1.1 延迟时间τ的确定

相空间重构对延迟时间τ并不十分敏感,τ只要选在合适的范围内即可,但这不意味着τ可以随意选取,若延迟时间τ过小,吸引子不能完全打开,相空间中的相点会积压在同一位置,反之会使各变量相互独立,丢失部分特性。目前计算延迟时间的方法有很多,如自相关函数法、互信息法等。本文采用自相关函数法[18]进行计算,其优点是算法简单计算迅速。对采集到的时间序列,计算其自相关函数:

(9)

取自相关函数第一次降到C(0)/e时对应的τ值作为最优延迟时间,此时相空间中的吸引子完全打开,能展现较多的细节。

2.1.2 嵌入维数m的确定

若嵌入维数m较小,吸引子无法完全展开,反之会使计算量剧增,增加误差。目前确定嵌入维数的方法主要有虚假近邻法(FNN)、CAO法等。本文结合时间序列特点采用FNN方法[19]计算嵌入维数,其优点是计算简便且准确度较高。FNN方法通过计算d维和d+1维空间的欧式距离的平方来确定临近点数量。其确定嵌入维的基本思想是:当重构相空间的维数m小于嵌入维de时,原本不相邻的点由于投影到低维空间从而变成相邻的点;反之若m大于de,则本应相邻的点变成了不相邻的点。系统的吸引子从m小于de到m大于de这个过程中,近邻点数会突然变小,由此确定最佳嵌入维数m。虚假近邻法的计算公式为:

(10)

式中:x(t)为采集到的时间序列,xNN(t+mτ)为x(t)的最近邻点,Rm(t)为x(t)和xNN(t)之间的距离。当m增大到某一值时,fm(t)的值不会随着m增大而发生明显变化,此时对应的m值即为最佳嵌入维数m。

2.2 递归分析法

2.2.1 递归图

递归图[20]是用来分析时间序列周期性与混沌性的重要方法,其核心思想是将高维运动状态的轨迹映射到二维平面空间中,将动力学特征可视化。其数学表示为:

(11)

(12)

当Rij的值为1时,在递归图的(i,j)处呈现黑点,当Rij的值为0时,在(i,j)处呈现白点。且由定义可知,递归图是关于主对角线对称的。

建立递归图需要确定嵌入维数m、延迟时间τ和临界系数ε三个参数。τ和m分别根据上节中的方法进行计算。ε需要根据实际情况进行具体选择,如果ε过小,生成的递归图中黑点会非常密集,不容易体现细节特征;如果ε过大,递归图中会出现大量空白,从而失去很多细节。由文献[20]可知,阈值应使递归度在0.01～0.05范围内,因此本文取递归度在0.02附近时对应的阈值。

递归图的纹理结构包含着原始原动力学系统中相空间吸引子的变化趋势和规律。通过分析不同参数下递归图的纹理特征,可以揭示锥齿轮系统的动力学特性随不同参数的变化情况。

递归图的纹理结构主要包括4种:均匀结构、漂移结构、突变结构和周期结构。均匀结构往往出现在随机信号中;漂移结构往往是因为系统中某些参数在缓慢地变化,表现为向地渐变行为;突变结构是由系统参数剧烈变化引起的,表现为递归图中出现大量的黑块;周期结构往往在振荡系统中,表现为周期状的结构(如棋盘结构等)。

2.2.2 递归量化分析

递归图可以对系统的动力学特性做定性分析,而递归分量化析法[20]主要依据递归图中点的密度、对角线和垂直线段的特征对系统动力学特性进行定量分析。本文主要分析递归度、递归熵和捕获时间[25]。

(1)递归度(RR):指递归图中递归点的百分比(递归图中黑点与总点数的比值)。

(13)

(2)递归熵(ENTR):分析递归图中对角线长度分布的Shannon熵。

(14)

式中,p(l)为对角线分布概率密度。递归熵反映了系统的复杂程度,ENTR值越大,系统结构越复杂,原时间序列的混沌特性就越强。

(3)捕获时间(TT):反映时间序列在某种状态下变化的持续时间。

(15)

式中:P(v)指递归图中长度为v的水平直线的分布数量,vmin为定义的最短水平直线长度。

上述参数中,ENTR是基于对角线的量化指标,TT是基于垂直线段的量化指标。ENTR可以分析系统从混沌到周期或从周期到混沌的变化情况,而TT可以分析系统从混沌到混沌的变化情况。

3 结果分析

本文以啮合频率ωh和刚度系数ε为主要控制参数,重点研究这些参数的变化对齿轮系统动力学特性的影响。取锥齿轮的中点螺旋角βm=35°,法面压力角αn=20°,锥齿轮模型方程组中的各无量纲参数取值如下:

3.1 啮合频率对系统非线性动力学特征的影响

图2和图3分别为锥齿轮系统在不同ωh下的相图和递归图。从图2可以看出,随着ωh增大,系统在周期运动与混沌运动之间变化,吸引子表现出复杂的动力学特性,但始终在一定范围内运动。当ωh=0.30、0.55和0.70时,吸引子呈现椭圆形结构,系统均做周期运动,系统振幅范围分别为1.5～2.5、1.8～2.1和1.92～1.96;当ωh=0.90时,吸引子出现了扭曲的部分,此时系统仍然做周期运动,振幅范围为1.8～2.0;当ωh=1.3时,吸引子呈现多周期运动状态,振幅范围减小到了1.85～2之间;当ωh=1.65时,吸引子的轨线相互缠绕、交叉,虽然杂乱但仍在一定范围内运动,振幅范围为-2～5之间,此时系统做混沌运动,且混沌特征较明显;当ωh=1.85和2.0时,系统回归到准周期运动,吸引子振幅范围也减小到0.5～3.5和1.5～2.5之间。

(a) ωh=0.30 (b) ωh=0.55 (c) ωh=0.70 (d) ωh=0.90

(e) ωh=1.30 (f) ωh=1.65 (g) ωh=1.85 (h) ωh=2.00图2 不同ωh时的吸引子结构

(a) ωh=0.30 (b) ωh=0.55 (c) ωh=0.70 (d) ωh=0.90

(e) ωh=1.30 (f) ωh=1.65 (g) ωh=1.85 (h) ωh=2.00图3 不同ωh时的递归图

从递归图中可以看出,当啮合频率在0.3、0.55、0.9和1.3时,递归图中都只有平行于主对角线的45°斜直线,而斜线之间的间距不同,说明在这些工况下系统均作周期运动,但周期各不相同;ωh=0.7时,递归图中同时显示出与对角线平行的分段的短线和长线,说明此时系统做2周期运动;当ωh=1.65时,递归图中出现了有规律且复杂的图案,类似混沌系统的结构,此时系统运动状态较为复杂;当ωh=1.85和2.0时,递归图中只有对角线平行的长直线,且间距不同,说明在这两种工况下,系统均做准周期运动。

表1给出了递归量化分析结果,由表1可知,当递归度接近时,随着ωh的增大,ENTR和TT的值上下波动,当ωh=1.65时,ENTR和TT都达到了最大值,由此得出ωh=1.65系统运动状态较复杂。

表1 不同ωh时的递归量化分析结果

综上所述,锥齿轮系统动力学特性会随着ωh的变化而发生较大变化,当ωh=1.65时系统状态最复杂,此时系统做混沌运动。因此为了减小弧齿锥齿轮系统的冲击,保证系统运动的平稳性,工程中应尽量避免啮合频率在这一值附近。

3.2 刚度系数对系统非线性动力学特征的影响

图4和图5分别为ωh=1.65,ξm=0.05时,系统随ε变化得到的相图和递归图。由图4和图5可以看出,随着ε增大,相图和递归图都发生了较大的变化。图5中吸引子从单一且规则的圆形轨迹逐渐扭曲变形,最后又变得规则。当ε=0.5和0.6时,吸引子相互纠缠在一起,其余工况下均为单一的曲线,但运动范围大小不同,ε=1.0时,吸引子的振幅范围最小。

(a) ε=0.1 (b) ε=0.3 (c) ε=0.5

(d) ε=0.6 (e) ε=0.8 (f) ε=1.0图4 ξm=0.05时不同ε下的吸引子结构

(a) ε=0.1 (b) ε=0.3 (c) ε=0.5

(d) ε=0.6 (e) ε=0.8 (f) ε=1.0图5 ξm=0.05时不同ε下的递归图

随着ε的增大,递归图结构发生了较大的变化。当ε=0.1、0.3、0.8和1.0时,递归图中只有平行于主对角线的长斜线,结合对应的相图可得出,在这些工况下,系统作单周期或准周期运动;当ε=0.6时,递归图中除了有规则的平行于主对角线的结构外,还出现了垂直方向的直线,结合相图可知此时系统做多周期运动;ε=0.5时,递归图结构较特殊,出现了较短的平行于主对角线的斜直线,且不同短线间的距离不同,说明此时系统做混沌运动。综上所述,随着ε的变化,系统从周期运动到混沌运动最后又回到周期运动,ε=0.5时,系统运动状态最复杂。

图6和图7分别为ωh=1.65,ξm=0.07时,系统随ε变化得到的相图和递归图。从图6中能够看出,当ξm=0.07时,随着ε增大,系统动力学特性发生较大变化,相图中吸引子的轨线从单一且规则的圆形逐渐变得扭曲且互相重叠,最后又变得规则。图7中,ε=0.1、0.8和1.0时,递归图中只有与主对角线平行的斜直线,在这些工况下系统做稳定的周期运动;ε=0.3时,递归图中右上方有较多的黑点,而左下方是大量白色区域,反映了此时数据的非平稳性,此时系统含有偏移的趋势,处于不稳定的运动状态;ε=0.5时,递归图中呈现出规则的图案,且图中出现垂直线,结合相图分析可推测,此时系统为多周期运动;当ε=0.6时,递归图呈现散点和块状共存的结构。综上所述,系统在ε=0.1、0.8和1.0时的运动情况较简单,在其他几种工况下运动较为复杂。

(a) ε=0.1 (b) ε=0.3 (c) ε=0.5

(d) ε=0.6 (e) ε=0.8 (f) ε=1.0图6 ξm=0.07时不同ε下的吸引子结构

(a) ε=0.1 (b) ε=0.3 (c) ε=0.5

(d) ε=0.6 (e) ε=0.8 (f) ε=1.0图7 ξm=0.07时不同ε下的相图

表2和表3给出了弧齿锥齿轮系统在不同ε下的递归量化分析结果。可以看出,当递归度RR接近,ξm=0.05和0.07时,随着ε的变化,ENTR的值均上下波动,TT的值均先增大后减小。ξm=0.05时,ε为0.5和0.6时,ENTR的值较大,ε为0.5时,TT的值最大,可知此时系统的运动状态较复杂;ξm=0.07时,ε为0.6时,TT和ENTR的值最大,此时系统运动状态较复杂。

表2 ξm=0.05时不同ε的递归量化分析结果

表3 ξm=0.07时不同ε的递归量化分析结果

综上所述,系统的动力学特性随着ε的增大发生明显的变化,系统由稳定的周期运动经历混沌、周期、混沌运动后,最后又变为周期运动。因此当ωh和ξm值确定时,为了保证锥齿轮传动系统的平稳性,应该选择较高的ε。

3.3 阻尼系数对系统非线性动力学特征的影响

图8和图9为ωh=1.65,ε=0.2时不同ξm下的相图和递归图。图10和图11分别为ωh=1.65,ε=0.5时不同ξm下的相图和递归图。由图8可以发现,随着ξm的增大,吸引子结构从规则且单一的图案变得复杂,出现扭曲且交叠在一起的结构;当ξm进一步增大到0.1时,吸引子恢复单一且规则的圆形结构,相图中振幅范围从-3～5减小到1～3的范围。由图9可以发现,递归图中始终存在较长的斜直线,说明系统始终具有一定的确定性。当ξm=0.05和0.07时,递归图中出现了黑色的竖线,且ξm=0.07时的黑色区域更多。

(a) ξm=0.01 (b) ξm=0.05 (c) ξm=0.07 (d) ξm=0.10图8 ε=0.2时不同ξm下的吸引子结构

(a) ξm=0.01 (b) ξm=0.05 (c) ξm=0.07 (d) ξm=0.10图9 ε=0.2时不同ξm下的递归图

从图10可以看出,当ε为0.5时,随着ξm的增大,相图的轨迹先变复杂然后变规则,振幅范围从-3～5逐渐减小到1～3。当ξm=0.05时,相图中的吸引子扭曲折叠在一起,可知此时系统呈现较明显的混沌特征;当ξm=0.1时,相图中的轨迹变为单一且规则的圆形。从图11可以看出,当ξm=0.01和0.1时,递归图中只有长的斜直线,结合相图可知,此时系统做稳定的周期运动;当ξm=0.07时,递归图中呈现出规则的图案,且图中出现垂直线,结合相图分析可推测,此时系统为多周期运动;当ξm=0.05时,递归图结构复杂,表明此时系统运动情况复杂。

(a) ξm=0.01 (b) ξm=0.05 (c) ξm=0.07 (d) ξm=0.10图10 ε=0.5时不同ξm下的吸引子结构

(a) ξm=0.01 (b) ξm=0.05 (c) ξm=0.07 (d) ξm=0.10图11 ε=0.5时不同ξm下的递归图

从表4和表5可以看出,在递归度接近时,当ε=0.2和0.5时,TT和ENTR的值均随着ξm的增大,而呈现先增大后减小的规律,在ξm为0.05时TT值最大,ξm为0.07时ENTR的值最大,两者均在ξm为0.1时较小。

表4 ε=0.2时不同ξm的递归量化分析结果

表5 ε=0.5时不同ξm的递归量化分析结果

综合上述分析可以发现,在不同ε下,随着ξm的增大,系统的非线性动力学特征变化情况是类似的。ξm增大到0.1时,相图中均只有单一闭合曲线,递归图中始终只有与主对角线平行的直线,且ENTR和TT均为最小值,此时系统呈现周期运动特征。因此,当ωh和ε的取值确定时,ξm应取较大值以保证系统运行平稳。

综合3.1～3.3节分析可知,为了保证锥齿轮系统运行的稳定性,应使ξm、ε保持较高值,ωh避开某些特定区域。

4 结论

本文建立了弧齿锥齿轮系统的数学模型,并用龙格库塔法求出了不同参数下的时间序列,利用相空间重构、递归图和递归量化分析等非线性动力学方法,分析了系统在不同啮合频率、刚度系数和阻尼系数下的非线性动力学特性,为弧齿锥齿轮结构设计和运行参数优化提供参考。结论如下:

(1)相空间重构、递归图和递归量化分析方法可以定性、定量地分析弧齿锥齿轮系统的运动状态,从而识别系统的非线性动力学特征。弧齿锥齿轮系统的非线性动力学特征会随着啮合频率、刚度系数和阻尼系数的变化而发生较大变化。

(2)在刚度系数较高、阻尼系数较大且啮合频率远离特定啮合频率时,弧齿锥齿轮系统的相图与递归图结构较简单,振幅范围较小,此时系统会表现出较强的周期性,反之,相图与递归图的结构较复杂,系统的非平稳性、非线性特征会增强,表现出混沌特征。

(3)在刚度系数较高、阻尼系数较大且啮合频率远离特定啮合频率时,系统的递归熵和捕获时间值均较小,此时系统运行平稳性较强,这与定性分析结果一致。因此,为了减小冲击,抑制混沌现象导致的系统运行不平稳,扩大系统的稳定工作区域,弧齿锥齿轮系统在设计和运行时应使啮合频率远离特定的区域,选用较大的刚度系数和较大的阻尼系数。