向 玲, 高 楠, 唐 亮, 郭鵬飛
(華北電力大學 機械工程系,河北 保定 071003)
齒輪傳動系統是航空航天、高鐵運輸和風能利用等領域中重要的機械設備。近年來風力發電機作為清潔能源的發電設備已得到全球的推廣和利用,但由于此類設備常在苛刻的環境下工作,其故障率也逐步增加,且齒輪傳動故障位居前列。其中,行星齒輪傳動是風電傳動系統的關鍵部件之一,但由于其結構復雜、內外激勵因素較多,會存在一些非線性因素。因此,齒輪傳動系統的動力學行為及其機理的研究,對提高系統效率和延長其工作壽命有重要意義。多年來,眾多學者致力于齒輪系統的動力學研究。邱星輝等[1]從系統動力學建模及求解、動力學特性和優化設計等方面對目前風電行星齒輪傳動系統的研究進行了綜述;Kahrarman[2]建立了單級行星齒輪傳動系統的動力學模型并對其振動模態做了分析;Lin等[3]考慮多個自由度建立了行星齒輪傳動系統的模型并進一步求解分析了系統固有頻率和振型,討論了陀螺效應和時變嚙合剛度對其的影響。Sheng等[4]建立了行星齒輪傳動系統的平移-扭轉非線性動力學模型,以系統內外嚙合副齒側間隙的變化作為分岔參數,引入阻尼的影響,深入研究了不同齒側間隙變化下行星齒輪傳動系統的動力學特性。李晟等[5]以兩級行星齒輪傳動系統為分析模型,研究了內外激勵下系統的分岔與混沌特性。林騰蛟等[6]以多級齒輪傳動系統為研究對象,求解并分析了齒側間隙、轉速及負載等因素對整個系統的影響。Wang等[7]建立了風電齒輪傳動系統的純扭轉非線性動力學模型,討論了阻尼的變化對系統非線性行為的影響。田亞平等[8]以風力發電機混合輪系為建模對象,求解并分析了該系統的動力學行為。楊軍[9]對行星齒輪傳動系統固有振動特性進行了研究,分析了包含支承剛度和時變嚙合剛度在內的多個參數對固有頻率及振型的影響。佘凱等[10]考慮支承剛度和支承阻尼,建立了直齒輪副的非線性動力學模型,分析支承剛度變化下系統的動力學行為。李貴彥[11]以弧齒錐齒輪為研究對象,建立了8自由度的非線性動力學模型,以支承剛度和嚙合間隙為分岔參數對其進行分析,探討了參數變化下系統的狀態。Xiang等[12]建立了含摩擦的齒輪-軸承非線性動力學模型,分析了支承剛度和激勵頻率變化下系統的非線性動力學特性。
當前許多學者已經建立了不同齒輪傳動系統的動力學模型,并考慮了多個非線性因素,研究分析了參數的變化對系統的影響。本論文在此基礎上,以1.5 MW風電齒輪傳動系統為研究對象,建立一級行星齒輪傳動及兩級平行軸齒輪傳動系統的平移-扭轉非線性動力學模型,模型考慮了時變嚙合剛度、綜合嚙合誤差、齒側間隙和太陽輪支承剛度,進一步地討論和分析太陽輪支承剛度的變化對系統動力學特性的影響。分析結果為深入研究風電齒輪傳動系統的非線性動力學特性和故障機理等提供理論參考。
圖1所示為某風電齒輪傳動系統動力學模型示意圖,該模型由一級行星齒輪傳動和兩級平行軸齒輪傳動串聯構成,所有齒輪均為直齒圓柱齒輪。傳動系統輸入端為行星架,輸入扭矩為Tin,輸出端為高速端齒輪系,輸出扭矩為Tout。行星齒輪系包括太陽輪、行星架、行星輪和內齒圈,分別用s、c、pi(i=1,2,…,N)和r表示。平行軸齒輪傳動包含低速端齒輪傳動和高速端齒輪傳動,共包含4個齒輪,分別以g1、g2、g3和g4表示。為簡化計算,假設全部齒輪主體部分為剛體,齒牙為彈性體,行星齒輪在行星架均布且參數一致,傳動軸為剛性短軸且兩端軸承具有相同剛度和阻尼,忽略各輪系齒間摩擦力的影響。行星架(各齒輪)的扭轉位移表示為θm(m=s,pi,c,g1,g2,g3,g4),太陽輪徑向位移為xs、ys,建立8+N自由度動力學模型,其廣義坐標為即
q={θs,xs,ys,θc,θp1,…,θpi,θg1,θg2,θg3,θg4}
齒輪內部激勵包括時變嚙合剛度、嚙合阻尼和嚙合誤差,分別由kj、cj、ej(j=spi,rpi,g1g2,g3g4)表示,其中spi代表行星齒輪傳動內嚙合副,即太陽輪和各行星輪嚙合副、rpi代表行星齒輪傳動外嚙合副,即內齒圈和各行星輪嚙合副、g1g2和g3g4則分別代表平行軸齒輪傳動中的低速端齒輪嚙合副和高速端齒輪嚙合副。

圖1 風電齒輪傳動系統動力學模型
齒輪嚙合副剛度的變化是時變的,根據直齒輪嚙合副剛度的特點可近似的采用周期矩形波來表示各齒輪副間嚙合剛度。為進一步簡化計算,本文采用以基于嚙合頻率ωj為基頻的Fourier級數展開表示時變嚙合剛度,取一次諧波項,公式可表示為[13]
kj=kmj+kajsin(ωjt+φj)
(1)
式中,kmj、kaj分別為各嚙合副的平均嚙合剛度及剛度變化幅值,φj為各嚙合副嚙合剛度變化幅值的初相位。
為使方程能夠順利求解,引入齒輪副間的相對位移xj,公示表達為
xg1g2=rg1θg1+rg2θg2-eg1g2(t)
xg3g4=rg3θg3+rg4θg4-eg3g4(t)
xrpi=rpiθpi-rcθccosα-erpi(t)
(i=1,2,…,N)
(2)
式中,ej(t)為綜合嚙合誤差,是靜態傳遞誤差和齒輪制造誤差的總稱,同時也是系統內部激勵來源之一,公式以嚙頻周期性變化的正弦函數來表示[14]
ej=Ejsin(ωjt+ηej)
(3)
式中:Ej為各嚙合副綜合嚙合誤差幅值;ηej為各嚙合副綜合誤差初相位。
齒側間隙非線性函數[15]表示為
(4)
系統中各嚙合副動態嚙合力由彈性恢復力和阻尼力構成,可表示為
(5)
式中,Mj是各齒輪副間的等效質量。

(6)
最終得到該系統的相對位移坐標下的無量綱振動微分方程表達式如下所示
(7)
式中,
該方程采用4~5階變步長Runge-Kutta法進行求解,為了消除系統的瞬態響應,略去前1 500個周期,計算所得數值結果以分岔圖、最大Lypapunov指數圖(以下稱LLE圖)、時間歷程圖、FFT頻譜圖、相圖及Poincaré截面圖作為分析工具對系統分岔及混沌特性進行研究說明。以兆瓦級風機作為示例,取齒輪傳動系統基本參數Pin=1.5 MW,Λ=94,ξ=0.03,α=20°,N=3,齒輪基本幾何參數和計算參數如表1所示。


圖2 系統隨激勵頻率變化的分岔圖(ki=1)

行星齒輪系行星架c太陽輪s行星輪p內齒圈r分度圓半徑rm/m0.470.190.280.75轉動慣量Im/(kg·m2)85.35.7529.6平均嚙合剛度kmj/(N·m-1)1.31×10101.486×1010剛度變化幅值kaj/(N·m-1)4.96×1095.12×109壓力角/(°)20平行軸直齒輪系低速端g1低速端g2高速端g3高速端g4分度圓半徑rm/m0.530.110.370.09轉動慣量Im/(kg·m2)3060.5939.90.17平均嚙合剛度kmj/(N·m-1)3.76×1095.12×108剛度變化幅值kaj/(N·m-1)7.25×1081.33×108

圖3 系統隨激勵頻率變化的最大Lyapunov指數圖(ki=1)
結合兩圖分析可知,隨著激勵頻率的增大,系統存在豐富的非線性行為,期間出現多個跳躍點和分岔點,系統多以激變和倍分岔的途徑進入混沌運動。
當激勵頻率Ω在范圍0.001~0.413 5時,系統處于單周期或擬單周期運動,期間伴隨陣發性混沌運動,對應的LLE也出現突變現象。當Ω在0.413 5~0.439的范圍內,系統進入混沌,LLE大于0。隨后系統回歸至擬單周期運動,在Ω=0.46時,系統分岔為擬2周期運動。當Ω=0.464 9時,發生跳躍現象,然后系統經激變進入混沌運動。隨著激勵頻率的增大,系統逐漸變化為擬2周期,在Ω=0.521 5時,發生又一次跳躍現象,系統進入短暫的單周期運動后發生了再一次跳躍現象后重新進入混沌。當Ω=0.589時,系統經激變進入單周期,經倍分岔后系統進入混沌運動,且混沌范圍明顯增大,對應的LLE也呈現增大的趨勢,僅在激勵頻率Ω=0.872 5時,系統由混沌激變為2周期,對應的LLE也出現突降,隨后經短暫的倍分岔后回歸混沌。
圖4所示為系統在支承剛度比ki=0.5時,系統隨激勵頻率變化的的分岔圖。對比圖2可知,當支承剛度減少時,系統總體運動情況不變,仍然以擬周期運動和混沌運動為主,但混沌區域明顯增強,如激勵頻率Ω在0.445~0.505的區域。而且,在Ω=0.8附近的周期窗口也明顯減小。

圖4 系統隨激勵頻率變化的分岔圖(ki=0.5)
圖5是支承剛度比ki=3時,系統隨激勵頻率變化的分岔圖,與ki=1和ki=0.5不同的是,系統在激勵頻率處于0.302 5~0.325的區間內,出現了混沌運動。對比圖2和圖4可知,系統的混沌運動有所抑制,出現了多個多周期窗口,如范圍在0.758 5~0.784和0.902 5~0.975 5等。在當激勵頻率Ω在0.758 5~0.784的范圍內,系統由混沌經倒分岔進入2周期運動,后又經激變進入混沌。當激勵頻率Ω在0.902 5~0.975 5的范圍內,系統變化狀態為混沌→8周期運動→6周期運動→2周期運動→4周期運動→8周期運動→混沌。

圖5 系統隨激勵頻率變化的分岔圖(ki=3)
以太陽輪支承剛度ks作為分岔參數,取剛度比ki在0.001~10區間,激勵頻率Ω=1時進行分析,得到全局分岔圖(含局部細化圖)及LLE圖如圖6和圖7所示??芍?,ki在0.001~0.081的范圍內,系統處于2周期運動,LLE小于0,其中ki在0.097~4.977區間,系統經激變進入混沌運動,期間出現多個倒分岔現象或短暫的多周期窗口,如在1.281~1.313、1.953~2.241和2.417~2.449的范圍等,可見對應的LLE出現突降現象。隨后,系統由混沌進入2周期運動,經分岔為4周期運動后又回歸至2周期,當ki=8.193時,系統再次發生分岔為4周期運動,直到ki=9.569,系統經激變重新回到2周期運動并保持不變。由圖7可知,在期間發生了多次陣發性混沌運動,其對應的LLE指數也是大于0的。

圖6 系統隨支承剛度比ki變化的分岔圖及其細化圖
圖8~圖12表示了系統經倒分岔由混沌進入2周期運動的變化過程,分別以系統在不同支承剛度下的時間歷程圖、頻譜圖、相圖和Poincaré截面圖做以說明。當ki=1.665時,系統處于混沌運動,如圖8所示。在相圖中相軌跡是無規則的運動,Poincaré截面為無規則的點集,時間歷程圖無明顯規律可循,對應的頻譜圖出現Ω/2分頻成分并存在明顯的連續邊帶。圖9所示為系統在ki=1.905時所呈現的混沌運動,對比圖8,相圖無規則軌跡逐漸聚攏且以內八軌跡為主,Poincaré截面有兩個點集,頻譜圖上出現Ωn/2分頻成分,n為正整數,不同的是Ω/2分頻高于基頻并存在些許連續邊帶。當ki=1.985時,如圖10所示,系統為8周期運動,頻譜圖上出現Ωn/8分頻成分,相圖軌跡為8個封閉曲線,Poincaré截面圖為8個點。當ki=2.081時,系統處于4周期運動,如圖11所示。在頻譜圖中可見其Ωn/4的分頻成分,其中Ω/2分頻幅值最高,相圖為4個封閉軌跡,Poincaré截面為4個點。當ki=2.161時,系統為2周期運動,如圖12所示。同樣,在頻譜圖中存在Ωn/2分頻成分,Poincaré截面為2個點,相圖軌跡呈現為2個交叉的封閉曲線。

圖8ki=1.665時系統的時間歷程圖、頻譜圖、相圖及Poincaré截面圖
Fig.8 The time series, FFT spectrum, phase diagram and Poincaré map whenki=1.665


圖9ki=1.905時系統的時間歷程圖、頻譜圖、相圖及Poincaré截面圖
Fig.9 The time series, FFT spectrum, phase diagram and Poincaré map whenki=1.905


圖10ki=1.985時系統的時間歷程圖、頻譜圖、相圖及Poincaré截面圖
Fig.10 The time series, FFT spectrum, phase diagram and Poincaré map whenki=1.985


圖11ki=2.081時系統的時間歷程圖、頻譜圖、相圖及Poincaré截面圖
Fig.11 The time series, FFT spectrum, phase diagram and Poincaré map whenki=2.081
綜上,支承剛度會對系統產生影響,同時表現出豐富的非線性動力學行為。支承剛度的提高會使系統出現由周期進入混沌的現象,但繼續增大則仍會回歸至周期狀態。因此在設計或作業時,應避開支承剛度的分岔點和混沌區域。


圖12ki=2.161時系統的時間歷程圖、頻譜圖、相圖及Poincaré截面圖
Fig.12 The time series, FFT spectrum, phase diagram and Poincaré map whenki=2.161
在風電齒輪傳動系統中,齒側間隙等的非線性因素可能會引起整體系統的不穩定性,而持續工作狀態中的支承剛度也會出現變化,從而導致系統出現分岔行為和混沌運動,破壞系統的穩定性,而混沌運動是系統不斷地由某種軌道突跳到另一個軌道上去的無規則運動,表示了系統不可預測的行為和軌道的永不重復性,會導致系統壽命減少、增加疲勞作業時間和提高系統磨損的可能性,進而產生噪聲,嚴重時則會致損。在設計中,應避開混沌狀態的范圍,采取安全、符合要求的作業范圍。
通過全局分岔圖、時間歷程圖、FFT頻譜圖、相圖及Poincaré截面圖可以定性對系統進行分析,利用全局或單點最大Lyapunov指數圖可以定量的對系統進行分析,從而可以更為精確的確定系統的運動狀態,保證機器的正常運行。
(1) 考慮齒輪非線性因素,如時變嚙合剛度、綜合嚙合誤差和齒側間隙,建立了風電齒輪傳動系統的非線性動力學的平移-扭轉模型,引入支承剛度作為分岔參數,分析了系統在不同支承剛度下隨激勵頻率變化的動力學特性。支承剛度的增大會對系統隨激勵頻率變化下的運動狀態產生影響,使進入混沌延遲并出現周期運動窗口。
(2) 利用多種分析工具對系統響應結果進行分析,在時變嚙合剛度、綜合嚙合誤差和齒側間隙強非線性因素的作用下,系統在太陽輪支承剛度的變化下會表現出豐富的分岔特性及混沌現象,如周期運動、擬周期運動和混沌運動,進入混沌運動的途徑以激變和倍周期分岔途徑為主。本文分析結果可對風電齒輪傳動系統的設計、安裝及減振方面提供理論參考。