直齒錐齒輪分岔脫嚙特性參數解域界結構

田亞平, 楊江輝, 王瑞邦

(蘭州交通大學 機電工程學院, 蘭州 730070)

0 引 言

錐齒輪傳動廣泛應用于機車、航空、機床等行業傳動裝置中．含間隙和時變因素的錐齒輪傳動系統的非線性振動對系統運動的穩定性提出了嚴峻挑戰,成為學者們研究的熱點．王三民等[1]建立了正交錐齒輪系統動力學模型并分析了系統分岔和幅值跳躍特性．黃康等[2]考慮隨機因素探索了間隙的合理范圍．Wang等[3]考慮含齒面混合潤滑研究了表面粗糙度對齒面沖擊的影響．Hua等[4]建立了軸承-錐齒輪轉子有限元模型,分析了軸承支承剛度對動態嚙合力的影響．Cao等[5]分析了不同潤滑狀態下錐齒輪傳動系統的非線性動態特性．李飛等[6]研究了摩擦和重合度對系統振動幅值和嚙合力的影響．

小波法[7]、數值法[8]是求解非線性系統分岔特性的常用方法,但在全局分岔特性轉遷特性分析上適用性較差．胞映射理論是研究非線性動態特性轉遷的常見方法．在參數域平面內,Gou等[9]采用胞映射研究了齒輪副模型的周期吸引子分形結構,Liu等[10]通過對齒輪解域界分析實現了混沌控制,林何等[11]分析了弧齒錐齒輪傳動系統中周期運動的解域界結構．因偽不動點周期求解及追蹤CPNF法求解效率高,在含間隙的行星輪系[12]和單級齒輪系統[13]的周期求解中獲得了成功應用．田亞平等[14]用CPNF法研究了單級齒輪系統的動態特性解域界結構．多參耦合下的解域界結構能有效揭示其動態特性間的耦合及轉遷規律．目前,針對含間隙直齒錐齒輪動態特性參數解域界結構的研究還鮮有報道．

據統計,間隙和時變參數是導致動車齒輪傳遞裝置(圖1)故障和車輛振動超標的主要因素．本文以該錐齒輪傳動系統為研究對象,建立了7自由度錐齒輪系統動力學模型,采用CPNF法求解周期、脫嚙、齒面沖擊、動載特性的轉遷規律,探索參數組合下系統的周期分岔與齒面沖擊及脫嚙、動載間的耦合關系,為錐齒輪結構設計參數選擇提供理論依據．

1 動力學模型

含彈性支承的直齒錐齒輪傳動系統動力學模型如圖2所示[1]．兩齒輪軸線正交于原點O,建立坐標系∑:{Oxyz}．兩輪支承等效在齒寬中點O1,O2,沿坐標軸的支承剛度和阻尼為kij,cij(i=x,y,z;j=1,2)．ch,kh(t)和bh分別為齒輪副嚙合的阻尼、時變剛度和齒側間隙．設主動輪1受到的驅動力矩Tp由常量Tpm和變量Tpv構成;從動輪2受到穩定的阻抗力矩Tgm．忽略齒輪的擺動與彎曲,兩個剛性錐盤的振動是分別以O1,O2為中心,沿坐標軸的移動和軸線的轉動,即8個自由度{X1,Y1,Z1,X2,Y2,Z2,θ1,θ2}．

考慮振動和誤差,兩齒輪沿嚙合線的動態綜合傳遞誤差Λn為

Λn=(X1-X2)a1-(Y1-Y2)a2-(Z1+r1θ1-Z2-r2θ2)a3-en(t),

(1)

式中,a1=cosδ1sinαn,a2=cosδ1cosαn,a3=cosαn,δ1為主動輪節錐角,αn為法面壓力角,r1,r2為兩輪的節圓半徑．齒輪副的靜態綜合誤差en(t)展開成級數形式為

(2)

式中,Ωh,Al和Φl為嚙合頻率、l階諧波幅值和初相位．

齒輪副的嚙合力Fn及其沿坐標軸的分力Fx,Fy和Fz分別為

(3)

用文獻[15]計算kh(t)并用級數展開為

(4)

式中,km,kkl和Φkl為嚙合剛度的均值、l階諧波幅值和初相位．結合動車齒輪箱工程計算精度要求,本文選取一階級數進行數值仿真計算．

間隙函數f(Λn)為

(5)

圖2 動力學模型Fig. 2 The nonlinear dynamic model for a spiral bevel gear set

根據Newton第二定律,圖2所示的錐齒輪傳動系統振動方程為

(6)

式中,m1,m2,J1,J2為兩輪的質量和轉動慣量．

用沿嚙合線的動態綜合傳遞誤差Λn作為新自由度,將式(6)中兩齒輪扭振方程合并得

(7)

將上述方程整理并量綱一化處理得

(8)

式中

xj=Xj/bn,yj=Yj/bn,zj=Zj/bn,λ=Λn/bn,Ωn=(km/me)1/2,

Ωij=(kij/mj)1/2,ξij=cij/(2mjΩn),kij= (Ωj/Ωn)2,i=x,y,z,j=1,2,

(9)

bn為標稱間隙．

2 系統動態特性仿真指標

2.1 脫嚙占空比和齒背嚙合比

輪齒脫嚙和齒背接觸行為是引起齒面拍擊和附加動載荷的主要因素．由文獻[16]知,齒輪嚙合的沖擊和脫嚙(I=0,1,2)狀態可由運動周期中的最小嚙合綜合誤差λmin與齒側間隙b判斷,其判斷式為

(10)

脫嚙占空比(non-meshing duty coefficient,δNMDC)和齒背嚙合比(back-meshing duty coefficient,δBMDC)是量化系統脫嚙沖擊嚴重性的主要指標,設在一個運動周期T內脫嚙時間為t1,齒背嚙合時間為t2,則δNMDC和δBMDC為

(11)

2.2 動載系數

齒輪強度設計中考慮齒輪振動導致動態嚙合力瞬態增大的因數稱為動載系數(dynamic load coefficient,δDLC)．受支承剛度、箱體結構和重合度等因素影響動載系數計算十分復雜,為突出間隙和時變剛度對動載系數的影響,定義運動周期內最大嚙合力|fn(τ)|max與平均值fpm之比的平方根為動載系數:

(12)

3 參數解域結構求解

參數解域結構求解是基于胞映射原理將參數平面進行胞元離散,在胞元內對周期運動采用CPNF法延續追蹤判穩獲得其周期運動的參數域界結構．在追蹤過程中伴隨有脫嚙、沖擊和動載特性的參數解域結構求解．其基本思想及實施過程如下．

3.1 參數域胞元化

對多參γi控制的非自治動力學系統可描述為

(13)

式中,X為n維狀態向量,τ為量綱一時間,γi∈RK為系統的k個控制參數,f為映射法則．

為分析激勵參數γi與響應間的關聯機制,在[xl,xu]×[yl,yu]二維參數域Θ內對參數離散化為n×m個胞元,參數胞為cij(i=1,2,…,n;j=1,2,…,m)．其中,xl,yl和xu,yu分別表示x,y方向的起點和終點,hi,hj為離散胞在x,y方向的胞元尺度．x,y方向的胞尺度hi,hj為

(14)

在參數域Θ內,從xl,yl開始得到x,y方向上參數的遞推關系xi,yj為

(15)

則域Θ內離散后的各胞元(xi,yj)賦值為

(16)

3.2 動力學行為延拓追蹤步驟

參數域Θ內各胞元的動力學行為采用CPNF法追蹤求解．首先沿x方向對初始胞元(x0,yj)用偽不動點PNF法[10-11]進行系統運動周期求解并判穩; 然后沿x方向按步長hi進行延拓追蹤、用Floquet乘子|λ|max判斷周期運動的穩定性．當|λ|max>1時周期運動失穩,在該胞元(xi,yj)處終止周期追蹤判斷分岔類型,并用PNF法重新求解新周期,若求解周期數超過設定值時認定該胞元為混沌或擬周期狀態(可用Floquet乘子和Lyapunov指數判定),在(xi+1,yj)胞元處重新求解周期數并追蹤直至xu結束;當完成x方向追蹤后,沿y方向遞進一層繼續追蹤求解,直至遍歷全域結束．

3.3 胞元周期解求解過程

CPNF法周期求解判穩的關鍵是求方程組(17)的穩定解

(17)

式中,Φ=DP(X)為系統的轉遷矩陣[12], 其狀態方程的Jacobi矩陣?f(X,τ)/?X在λ≠±b的分段光滑區內為

(18)

式中

A1=[-a1[kx12ξx1],a2[ky12ξy1],a3[kz12ξz1],

a4[kx22ξx2],-a5[ky22ξy2],-a3[kz22ξz2]],

(19)

因在λ=±b處Jacobi矩陣不存在,故結合計算精度選取合理步長的有限差分法替代為

(20)

方程(17)以(Xk,I)為初值積分一個周期T獲得系統迭代點Xk+1和轉遷矩陣DP(Xk)．令Xk=Xk+1作為下次迭代的初始值,迭代到滿足精度要求的不動點XF,就是式(13)的m倍周期運動X(τ)的不動點．找到穩定周期解后,引入分岔參數Ω,采用延續追蹤法進行周期追蹤．設已求得Ω=Ωk(k=0, 1, 2, …)時,周期T的不動點為Xk,用Euler積分法得Ω=Ωk+1點對應不動點Xk+1的初始不動點Xk+1,0的預測公式[14]:

(21)

式中,ΔΩ為追蹤步長,GX(Ωk,Xk)=I-DP(Ωk,Xk),GΩ(Ωk,Xk)=PΩ(Ωk,Xk)．PΩ(Ωk,Xk)由打靶法建立的微分方程組(22)求得

(22)

式中

fΩ(X,Ω,τ)=feΩ{2cos(Ωτ)+Ωsin(Ωτ)}[01×131]T+

αsin(Ωτ)[0 -1 0 1 0 1 0 1 0 -1 0 -1 0 -4]T/4．

(23)

以(Xk,Ωk)為初值積分一個映射周期T獲得以(Xk,PΩ(Ωk,Xk))為追蹤參數下的新不動點．

4 系統動態特性綜合分析與解域界

不失一般性,取量綱一仿真參數:

ξi1=ξi2=0.01 (i=x,y,z),ξh1=ξh2=0.012 5,ξh=0.05,

ki1=ki2=1.0,kh1=kh2= 0.5,kh=1+αcos(Ωτ),

fpm=0.5,fav=0.0,fe= 0.2,b=1.0．

在參數解域界平面內,齒面沖擊/運動周期(I/P)類型用特定的顏色表示,其分岔曲線為解域界．齒輪系統常常穩定運行于短周期狀態下,為突出低周期(P<64)分岔特性,在I/P分岔圖中對周期數P≥65的長周期、擬周期(quasi-P)和混沌(chaos)運動均用“I/N”表示．系統的倍化(PD bif)、鞍結(NS bif)、檫切[17](G bif)分岔定義參見文獻[14],不再贅述．當系統從I/P向I/(3P)狀態轉遷時,系統發生了齒輪系統典型的周期3分岔(3T bif),依次類推．

時變嚙合剛度系數α表征齒輪嚙合重合度,其值表征單雙對齒嚙合不均勻性程度．頻率比Ω表征系統轉速,其值表征系統轉速高低．故構建Ω×α∈[1.2, 1.7]×[0.05, 0.55]參數平面(圖3),在平面內均勻地選取501×501個胞元構建I/P及相應的δNMDC,δBMDC,δDLC解域界結構．

I/P解域界結構(圖3(a))表明:在參數平面內,以周期1(P=1,下同)運動為主,存在P=2,3,9,10,18,32和quasi-P、chaos等運動狀態;3種齒面沖擊共存,以單邊沖擊為主．在α<0.15區域內,隨Ω遞增系統發生了SN bif和3T bif,使系統處于0/1,1/1,1/3運動狀態．在0.15<α<0.3區域內,系統通過PD bif、3T bif、HP bif等分岔方式進入混沌運動．在0.3<α<0.55區域內,隨Ω遞增系統發生了G bif、3T bif、3T chaos和CIC bif 等分岔．在Ω=1.6附近表現出復雜的動力學行為轉遷,即在該胞元內動力學行為極其不穩定．沿α方向在Ω<1.25區間內,系統以P=1運動的3種沖擊轉遷為主．在1.25<Ω<1.5區間,系統的動力學行為較為穩定．當Ω>1.5時,系統通過HP bif、3T bif和CIC bif的形式使系處于周期-混沌-周期運動狀態．參數解域界結構表明齒輪嚙合重合度是影響齒輪沖擊的主要因素,隨α遞增齒面脫嚙、齒背接觸的概率增大．在系統的共振頻率附近多周期、混沌現象加劇,并出現了齒輪傳動系統典型的周期3分岔行為．

系統I/P解域界結構對應的齒面脫嚙、齒背接觸和動載特性如圖3(b)—3(d)所示．圖3(b)表明,齒輪脫嚙程度隨周期分岔和齒面沖擊的轉遷而階躍突變,在1/1解域內其脫嚙程度并不一致,在G bif分岔線附近δNMDC達到極值0.55左右．隨α增大其脫嚙程度加劇,而在1/N和2/N解域內δNMDC值并未達到極值,即多周期或混沌運動下齒面脫嚙未必最嚴重．圖3(c)表明,在雙邊沖擊的2/1域界內其齒背接觸最為嚴重達到極值0.38左右,在2/3、2/N域界內其齒背接觸狀態隨α遞增而加劇但沒有達到極值,即齒輪嚙合的雙邊沖擊程度與系統混沌現象并無必然聯系．圖3(d)表明,系統的δDLC隨齒面沖擊和周期運動分岔出現了突變,在同一解域內其動載系數隨α遞增而增大,在無沖擊的0/1解域內其值較小、在雙邊沖擊的2/1、2/N域界內達到極值2.8．

(a) I/P(b) δNMDC

(c) δBMDC(d) δDLC圖3 Ω×α平面動態特性解域界結構Fig. 3 Dynamic properties solution domain boundary structures in the Ω×α 2-parameter plane

為驗證CPNF法求解的正確性,沿α方向以0.05等距選取截面,用Runge-Kutta法對Ω參數分岔進行數值計算獲得三維I/P分岔圖(圖4)．

圖4 Ω×α雙參分岔圖Fig. 4 The bifurcation diagram in the Ω×α 2-parameter plane

圖4中,黑、紅、藍色映射點分別對應無沖擊、單邊沖擊和雙邊沖擊,其I/P分岔過程與圖3(a)相匹配驗證了解界域結構的正確性．因CPNF法對周期運動求解追蹤過程不用計算瞬態,故較Runge-Kutta法求解效率高．

為更清晰地展示周期分岔與齒面脫嚙的關聯關系,取α=0.20,0.40的分岔圖如圖5所示．圖5(a)中的系統通過SN bif、3T bif、HP bif和激變(CIC bif)分岔方式使系統處于0/1, 1/1, 2/3, 1/3, quasi-P,2/chaos,1/chaos,1/1運動．圖5(b)中系統通過SN bif、3T bif、PD bif和C bif分岔使系統處于2/1,1/1,2/3,2/chaos運動狀態．HP bif和SN bif可通過CPNF追蹤中Floquet乘子|λ|max來判定,而G bif與SN bif是通過齒面沖擊類型和Floquet乘子共同判別,若穿越λ=±b截面時Floquet乘子從實數軸方向穿越單位圓則發生了鞍結分岔,否則為擦切分岔．

取α=0.20截面(圖5(a))的部分Poincaré映射和相圖(圖6)驗證CPNF求解的分岔過程．Ω=1.531時系統為1/1運動(圖6(a))經3T分岔為2/3運動(圖6(b)),P=3轉遷為3個焦點(圖6(c))、Hopf圈(圖6(d))、扭曲鎖相(圖6(e))直至形成混沌運動(圖6(f))．即圖5驗證了通過周期3分岔和Hopf分岔使系統處于1/1,2/3,1/1,擬周期和混沌運動．

(a) α=0.20(b) α=0.40圖5 分岔圖Fig. 5 The bifurcation diagram via α

(a) Ω=1.531(b) Ω=1.535(c) Ω=1.563

(d) Ω=1.575(e) Ω=1.600(f) Ω=1.650圖6 Poincaré映射圖和相圖Fig. 6 Poincaré and phase maps

圖7為圖5(b)中α=0.40時I/P部分的Poincaré映射和相圖．Ω=1.523時系統處于0/1運動(圖7(a))經周期3分岔為2/3運動(圖7(b)),P=3再經周期3分岔為P=9(圖7(c))、P=18(圖7(d))、……、混沌運動(圖7(e))．圖7驗證了系統經周期3分岔進入混沌的過程．

取α=0.40截面(圖5(b))對應的齒面脫嚙比、齒背接觸比和動載系數分布圖(圖8)．齒面脫嚙在齒面沖擊轉遷截面發生了突變,在1/1周期區間其脫嚙程度隨Ω遞增而減弱(圖8(a))．齒背接觸在雙邊沖擊區域隨Ω遞增而加劇,在陣發性混沌區域存在突變現象(圖8(b))．動載系數隨齒面沖擊和周期分岔產生了突變,隨齒面沖擊加劇其動載系數俱增(圖8(c)),在單邊沖擊區其δDLC隨Ω遞增而遞減,在雙邊沖擊區隨Ω遞增而遞增．即齒面沖擊和周期分岔是影響齒面脫嚙、齒背接觸及動載特性的主要因素．

圖8 脫嚙比、動載系數分布圖Fig. 8 The δNMDC, δBMDC and δDLC distribution diagrams

時變嚙合剛度和頻率比平面內的I/P,δDLC,δBMDC和δNMDC解域界結構揭示了系統隨參數分岔、脫嚙、沖擊、動載荷的轉遷規律及其關聯關系．針對不同速度車型要求,依據解域界結構對轉速和重合度參數進行匹配優化．選取參數時盡量避開混沌、雙邊沖擊和動載較大的參數區域,盡可能在0/1區域選擇合理的重合度或轉速參數能有效地降低嚙合沖擊、脫嚙和動載對齒輪箱的影響,可延長齒輪箱的服役周期、降低設備維護費用．

5 結 論

本文針對7自由度直齒錐齒輪傳動系統的周期分岔、齒面沖擊、脫嚙和動載特性問題,采用胞映射理論構建了重合度和轉速激勵參數平面,用CPNF法進行了數值求解獲得了I/P,δNMDC,δBMDC和δDLC解界域結構,得到了如下結論:

1) 參數平面內基于胞映射理論的CPNF法解域界結構求解方法是一種高效的數值求解方法．

2) 含間隙的錐齒輪系統存在鞍結、倍化、Hopf、激變、擦切和周期3等分岔現象,隨時變嚙合剛度系數增加其齒面沖擊加劇,在嚙頻附近系統的周期運動和齒面沖擊最為嚴重．

3) 系統分岔和齒面沖擊轉遷導致其齒面脫嚙和齒背接觸及動載行為產生了突變,在各I/P區域內其值也隨參數波動．隨時變嚙合剛度系數增大其動態特性加劇．在1/1和2/N區域其脫嚙最為嚴重,在2/1區域內齒背接觸最為嚴重,其動載系數最大．