考慮時變嚙合剛度的風機傳動鏈動態性能*

劉 宏，趙榮珍，鄭玉巧

(1.蘭州理工大學數字制造技術與應用省部共建教育部重點實驗室 蘭州，730050)(2.蘭州理工大學機電工程學院 蘭州，730050)

?

考慮時變嚙合剛度的風機傳動鏈動態性能*

劉 宏1,2，趙榮珍1,2，鄭玉巧1,2

(1.蘭州理工大學數字制造技術與應用省部共建教育部重點實驗室 蘭州，730050)(2.蘭州理工大學機電工程學院 蘭州，730050)

針對大型風力發電機組齒輪傳動鏈動態剛度引起的機組結構振動問題，綜合輪齒彎曲變形、齒根過度圓角處的基體變形和接觸變形等因素，建立齒輪時變嚙合剛度的量化分析模型，并與有限元動態嚙合模型對比驗證理論模型的正確性。在此基礎上考慮齒輪時變嚙合剛度和軸扭轉剛度推導1.5 MW風力機傳動鏈的動態總剛度，用于分析傳動鏈在動態剛度下固有特性變化規律及傳動鏈臨界轉速對動態剛度參數的敏感性，量化顯示動態剛度幅值變化引起的臨界轉速波動。研究表明,齒輪時變嚙合剛度的波動會引起傳動鏈臨界轉速的不穩定，增大時變剛度幅值會引起轉子系統臨界轉速的升高，但總體上嚙合剛度波動對臨界轉速的影響處于非敏感區。本研究對揭示風力機齒輪傳動鏈的內部剛度激勵機理和實現系統動態性能優化設計提供理論依據。

風力發電機； 傳動鏈； 嚙合剛度； 臨界轉速； 參數敏感性

引 言

增速箱是風力機組的關鍵部件之一，常工作于重載、沖擊和變載荷等復雜工況下，且機組高空架設、維修較困難等問題對其可靠性提出了更高的要求[1]。傳動齒輪對的嚙合剛度是傳動鏈內部重要的激勵源，掌握其運行機理對于研究傳動鏈動力學特性和可靠性設計極其重要，因而圍繞齒輪嚙合剛度的風力機傳動鏈動力學分析和設計的研究受到了廣泛關注[2-3]。

目前，國內外針對風力發電機齒輪傳動鏈的研究尚處于起步階段[4]，相關研究主要是將齒輪嚙合剛度簡化為定值模型展開，而傳動鏈的動態響應受齒輪時變嚙合剛度的影響[5]，特別在高速重載情況下更為明顯。因此，為保證分析的精確性需要考慮時變嚙合剛度的傳動鏈模型。在這一問題的研究中有限元方法[6-7]應用較為廣泛，但該方法需要模型網格細化工作并占用較多的計算時間，且風力機齒輪傳動鏈結構復雜，內部嚙合齒輪對較多，難于細化模型。影響振動的因素，特別是振動與其他運行參數之間關系非常復雜。然而，大型風力發電機組通過實測的方法難度大、成本高，基于理論分析的方法考慮齒輪時變嚙合剛度建立風力機傳動鏈行星輪系部分的動力學模型[8]，通過仿真計算，能夠獲得齒輪動態嚙合力，為時變嚙合剛度的理論驗證提供新的思路。同時，一些學者將時變嚙合剛度計入傳動鏈模型中，分析系統的動態響應特性[9-10]和齒輪傳動鏈的振動水平，探討傳動鏈的優化設計方法，加深了對風力機傳動鏈的認識。

基于上述情況，本研究以風力機傳動鏈中單對嚙合齒輪為研究對象，對建立齒輪時變嚙合剛度的量化分析模型問題進行研究，采用有限元動態嚙合模型對比驗證，在此基礎上細化增速齒輪箱，推導1.5 MW風力機齒輪傳動鏈的動態總剛度，量化分析傳動鏈在動態剛度下的固有特性及傳動鏈臨界轉速對動態剛度參數的敏感性。

1 輪齒嚙合剛度的量化分析模型

1.1 彎曲變形

齒輪嚙合剛度反映了齒輪接觸過程中接觸齒輪對數目的周期性變化及其嚙合線從齒根到齒頂的不均勻變化過程。它是一個周期性變化量，包括齒形的幾何特征、接觸點位置、輪齒撓度及齒形誤差等多方面因素。依據輪齒的彎曲、齒根過度圓角處的基體變形和接觸變形計算單齒的剛度和齒輪對的嚙合剛度。本研究將輪齒看作如圖1所示的不均勻懸臂梁，齒根部分懸臂固定在基體上，梁的有效長度為Le，沿嚙合線方向將輪齒劃分為n個微段，取其中一微段ei，假設嚙合齒輪對之間的嚙合力為F，由嚙合原理可知，F的大小和方向隨著嚙合點的移動而變化。以δb表示輪齒的彎曲變形，則輪齒的彎曲變形可表示為各微段變形的總和[11]

(1)

其中：α為壓力角；G為材料剪切模量。

(2)

(3)

其中：Ii為極慣性矩；Ai為輪齒接觸段面積；E為彈性模量；ν為泊松比。

由梁的彎曲剛度的計算方法可計算輪齒的彎曲剛度，表示為

kb=F/δb

(5)

圖1 輪齒嚙合模型

1.2 齒根過渡圓角處的基體變形

根據Muskshelishvili理論在彈性圓環中的應用，假設齒根圓部分受到恒定的隨時間變化的作用力，則該變形量δf[12]為

(6)

其中：W為齒寬。

L*，M*，P*和Q*可利用多項式函數式求得

(7)

其中：sf為單齒齒根對應齒輪基體的弧長；θf為sf/2所對應的齒輪圓心角；hfi為齒根圓半徑與基體內孔半徑的比值；系數Ai，Bi，Ci，Di，Ei和Fi的取值見表1[12]。

齒根過度圓角處的剛度系數kf可表示為

kf=F/δf

(8)

表1 多項式函數系數取值

Tab.1 The coefficient of polynomial function

系數L*M*P*Q*Ai×10-5-5.574060.111050.9520-6.2042Bi×10-3-1.9986-1.9986185.50009.0089Ci×10-4-2.3015-83.43100.05384.0964Di×10-34.7702-9.92560.05384.0964Ei0.02710.16240.2895-0.1472Fi6.80450.90860.92360.6904

1.3 接觸變形

根據已有文獻的研究結論，接觸齒輪對的赫茲接觸剛度可近似表達為式(9)所示的常數kh[13]。其數值大小與接觸線位置和滲透深度無關，僅由齒寬和齒輪材料的力學性能決定

kh=πEW/4(1-ν2)

(9)

1.4 齒輪對的嚙合剛度計算方法

綜合單齒彎曲變形、齒根過度圓角處基體變形和接觸變形，嚙合齒輪對的嚙合剛度ksp可表示為

ksp=1/(1/kbs+1/kbp+1/kfp+1/kfs+1/kh)

(10)

其中：kbs和kbp分別為太陽輪和行星輪的單齒彎曲剛度；kfp和kfs分別為太陽輪和行星輪齒根過度圓角處基體變形剛度。

2 時變嚙合剛度的計算方法設計

2.1 方法設計

一對齒輪嚙合傳動的區間是有限的，為了能夠連續傳動，必須保證在前一輪齒尚未脫離嚙合時，后一對輪齒能及時進入嚙合。單齒和雙齒嚙合的周期性變化引起齒輪嚙合剛度的周期性波動，從而造成了系統內部振動和噪聲，其嚙合區域的變化如圖2所示。從齒輪嚙合原理的角度，分析影響輪齒嚙合剛度變化的因素，主要由以下方面造成：a．在單齒嚙合周期內同等嚙合力作用下，因輪齒截面的不均勻性，使得輪廓齒面各點位移不同；b．整個齒輪嚙合過程中，單雙齒嚙合周期性變化使得輪齒位移和嚙合作用力也是變化的。

以傳動鏈中行星輪系傳動建立輪齒嚙合模型，主軸將動力傳遞給行星架，由行星架帶動3個行星輪驅動中心的太陽輪，假設太陽輪的運行速度為n1，且有Z1個齒，行星輪齒數為Z2，則單個齒的嚙合周期Te=60/n1Z1。該對嚙合齒輪的重合度εα為

εα=[Z1(tanαa1-tanα′)+Z2(tanαa2-tanα′)]/2π

(11)

其中：α′為嚙合角；αa1和αa2分別為齒輪1和2的齒頂圓壓力角。

圖2中Lc為單齒嚙合周期，根據太陽輪和行星輪參數計算得齒輪嚙合重合度。將輪齒接觸區域劃分為3段，在接觸的起始端(a端)有2對輪齒接觸嚙合(a～c段)，然后1對輪齒退出嚙合(c～e段)，隨下一對輪齒相繼嚙合，進入雙齒嚙合區域(e～g段)。考慮整個齒輪連續傳動重合度εα，實際動力學方程中需要的是齒輪綜合時變嚙合剛度，因此整個齒輪的時變內核剛度是有多對輪齒剛度綜合而成，對齒嚙合剛度在嚙合周期內疊加獲得整個齒輪的嚙合剛度。

圖2 輪齒接觸區域

2.2 應用驗證情況

以1.5 MW風力機齒輪箱內太陽輪及與之嚙合的行星輪為研究對象，其具體參數見表2。計算得太陽輪和行星輪嚙合重合度εα=1.51，雙齒嚙合區域占總嚙合區域的33.8%，單齒嚙合區域占總嚙合區域的32.4%。

表2 齒輪參數

在機組額定風速工作條件下由主軸輸入到太陽輪的扭矩為179.4 Nm，將扭矩向嚙合齒面轉化，假設3個行星輪承受相同的載荷條件，則施加到每個行星輪輪齒上的力為1.2 kN，由式(5),(7),(8)和(9)借助Matlab計算1對輪齒從進入嚙合到脫離嚙合1個周期內的剛度值，計算結果如圖3所示。

圖3 單齒嚙合剛度

從圖3可知，輪齒嚙合剛度的大小與嚙合位置密切相關。從齒根到齒頂嚙合的過程中，靠近齒根部分區域嚙合剛度迅速減小，越過齒輪分度圓后又緩慢增大，但齒頂附近的最大值小于齒根部分嚙合剛度的最大值。通過單對齒嚙合剛度在嚙合周期內疊加獲得整個齒輪的嚙合剛度曲線見圖4。

圖4 齒輪時變嚙合剛度

從圖4可知，時變嚙合剛度是一個周期性變化的量，其變化周期不再是單齒嚙合周期Lc，而是被重合度εα重新劃分為(2-εα)Lc。整個周期內有兩種變化趨勢，且單、雙齒變化之間為突變過程，嚙合剛度的最大值出現在雙齒嚙合區，最小值出現在單齒嚙合區。

3 輪齒嚙合剛度的有限元模型構造

采用三維實體建立行星輪系齒輪嚙合有限元模型，由于在有限元齒輪靜態接觸模型中，對兩個嚙合齒輪施加了全約束，實體單元沒有旋轉自由度，不能反映嚙合齒面上各點隨嚙合運動的位移變量。結合齒輪的工作要求，建立齒輪動態嚙合模型，該模型結合Ansys的參數化設計語言APDL建立，包括完整的齒輪結構。約束控制節點3個方向的平動自由度，保證齒輪的嚙合運動，向主動輪施加扭矩。為獲得一個齒嚙合范圍內嚙合齒面沿嚙合線方向各點的位移和嚙合力大小，在齒面嚙合線方向上均勻取10個節點作為1個節點集，將這些節點位移定義為變量，從而獲取各嚙合節點在單齒嚙合周期內沿x和y兩個方向的變化量δx和δy，則總位移為

δ=δxcosα+δysinα

(12)

節點集的變形計算結果如圖5所示，則單齒嚙合剛度計算如下

k=F/δ

(13)

圖5 嚙合齒面節點位移

圖6給出了采用理論分析模型和有限元仿真模型兩種方法求解得到的單齒動態嚙合剛度結果曲線。對比分析顯示，有限元動態嚙合模型在Ansys軟件中所模擬的單齒嚙合剛度值與理論計算的嚙合剛度值在輪齒進入嚙合和脫離嚙合的很小區間內誤差較大，其余嚙合區間上剛度值基本相同，且具有相同的變化趨勢，證明了理論分析模型在計算單對齒輪嚙合剛度上的有效性。

圖6 嚙合剛度對比分析

4 考慮時變嚙合剛度的風力機傳動鏈動態響應

4.1 傳動鏈動力學建模

將動剛度模型應用于1.5MW水平軸風力機傳動鏈(包括1級行星級和2級平行級齒輪傳動)中，行星級動力由行星架輸入，太陽輪輸出，內齒圈固定的結構形式，并用3個行星輪承擔載荷。其動力學模型如圖7所示。

圖7 傳動鏈動力學模型

模型中齒輪和軸假設成剛體，根據有關研究結論，齒輪嚙合振動的軸向、徑向振動中，軸向扭轉振動最為重要[12]。因此，只考慮系統的扭轉振動，應用拉格朗日方法建立系統的動力學方程

(14)

其中：q={qcqpqsq4q5q6q7}T為廣義坐標表示各自由度處的位移；F為外部力矢量；J為系統質量矩陣，J=diag[Jc3Jp+3mdJsJ4J5J6J7]；K(t)為由系統齒輪嚙合剛度和軸扭轉剛度組成的時變矩陣。

(15)

其中：ki(t)為齒輪時變嚙合剛度；kj為軸扭轉剛度(i=c,p,4,6;j=s,5)。

系統的剛度矩陣是一個時變矩陣，可表示為一個常數矩陣(平均值)k1和一個動態矩陣(波動矩陣)k1(t)的線性疊加，即

K(t)=k1+k1(t)

(16)

工程上，一般在討論結構固有特性時不計阻尼作用，由運動方程考慮系數矩陣k1可求得其無阻尼特征值，于是式(14)可變為

(17)

其中：φi為系統的陣型矩陣；ωi為系統自由振動的圓頻率。

按照機組系統的設計參數，齒輪動態嚙合剛度按前述方法依次進行參數化建模獲取，計算系統的固有頻率見表3。

表3 系統固有頻率

4.2 傳動鏈臨界轉速分析

傳動鏈兩級平行軸部分轉速不一致，為分析方便，把系統各部分向行星輪系輸出軸當量化處理，以太陽輪軸為等效構件，將平行軸的力和位移均用等效構件上的等效量來代替，當量化前后保證系統的動能和勢能守恒。將齒輪看成圓盤，圓盤間用嚙合剛度相連，傳動軸以扭轉剛度替代，再將各圓盤視為支點的回轉圓盤，通過傳動軸傳遞轉動和扭矩，則系統成為單軸系統。由于風載荷的隨機性，將風輪的發電轉速(11 ～ 22 r/min)折算到高速軸上，在Ansys中建立不同角速度下轉子的多載荷步模態分析模型，并計算得到系統的轉速頻率分布圖，如圖8所示。

由坎貝爾圖解得到轉子系統的1階臨界轉速為2 735.04 r/min，2階為3 958.9 r/min，3階為7 306.7 r/min，4階為10 530.6 r/min。額定工況下系統的額定轉速為1 500 r/min，計算可得工作轉速偏離其臨界轉速的裕度為82%，工作轉速相對于臨界轉速的裕度遠大于30%，因此裕度較大，滿足工作轉速偏離臨界轉速裕度的規定[14]。

圖8 傳動鏈轉速-頻率分布圖

4.3 參數敏感性分析

為了進一步研究齒輪時變嚙合剛度與傳動鏈的耦合作用，進行齒輪時變嚙合剛度參數對系統臨界轉速敏感性分析。考慮時變嚙合剛度從平均值7.74×108N/mm增大至動剛度的最大值9.86×108N/mm(增大21.5%)時，再次對系統分析模型進行求解，將變化前后所得系統臨界轉速進行對比，分析參數敏感性，如表4所示。

表4 轉子臨界轉速對比

對比分析發現，嚙合剛度的變化對系統臨界轉速的影響并不明顯，隨著嚙合剛度波動值提高21.5%，達到波動的最大值時，系統1階臨界轉速增大4.1%，即嚙合剛度稍有變化不會引起臨界轉速的很大變化，故表現出對該參數的不敏感。因此，系統有相對確定的臨界轉速，有利于系統工作頻率避開固有頻率而穩定地工作。

5 結 論

1) 將輪齒視為懸臂梁，綜合輪齒彎曲變形、齒根過度圓角處的變形和接觸變形建立了齒輪時變嚙合剛度的量化分析模型，并采用Ansys參數化設計語言APDL建立動態嚙合模型進行驗證，對比顯示兩種方法的計算結果，驗證了時變嚙合剛度理論建模方法的可行性。

2) 基于此，建立了1.5 MW水平軸齒輪傳動鏈的動力學模型，獲取了傳動鏈動態總剛度矩陣，1階臨界轉速為2 735.04 r/min，工作轉速相對于臨界轉速有充分的裕度。

3) 通過繪制不同嚙合剛度下轉子系統的坎貝爾圖，分析了系統臨界轉速與嚙合剛度之間的量化關系，結果顯示時變嚙合剛度的波動會引起系統臨界轉速變化，但總體上呈現系統臨界轉速對動態剛度不敏感。

[1] 陳濤,孫偉,張旭.基于灰色關聯度的風電齒輪箱傳動系統故障樹分析[J].太陽能學報,2012,33(10):1655-1660.

Chen Tao, Sun Wei, Zhang Xu. Reliability analysis of gearbox transmission for wind turbine based on fault tree model [J]. Acta Energiae Solaris Sinica, 2012,33(10):1655-1660. (in Chinese)

[2] 張義民,何永慧,朱麗莎,等.多平行軸齒輪耦合轉子系統的振動響應[J].振動、測試與診斷,2012,32(4):527-531.

Zhang Yimin, He Yonghui, Zhu Lisha, et al. Vibration response of multi-shaft rotor system with coupled gear mesh[J]. Journal of Vibration, Measurement & Diagnosis, 2012,32(4):527-531. (in Chinese)

[3] 蘇春,周小荃.基于有效年齡的風力機多部件維修優化[J].東南大學學報,2012,42(6):1100-1104.

Su Chun, Zhou Xiaoquan. Maintenance optimization for multi-component of wind turbine based on effective age[J]. Journal of Southeast University, 2012,42(6):1100-1104. (in Chinese)

[4] Peeters J L M, Vandepitte D, Sas P. Analysis of internal drive train dynamics in a wind turbine[J]. Wind Energy, 2006,9(1):141-161.

[5] Helsen J, Vanhollebeke F, Deconinck F, et al. Insights in wind turbine drive train dynamics gathered by validating advanced models on a newly developed 13.2MW dynamically controlled test-rig[J]. Mechatronics, 2011,21(4):737-752.

[6] 朱麗莎,張義民,馬輝,等.直齒輪耦合轉子系統的振動可靠性研究[J].振動、測試與診斷,2013,33(2):258-262.

Zhu Lisha, Zhang Yimin, Ma Hui, et al. Iterative learning control on double-valve parallel electro-hydraulic servo force control system [J]. Journal of Vibration, Measurement & Diagnosis,2013,33(2):258-262. (in Chinese)

[7] 于印鑫,袁惠群,梁明軒,等.基于顯式動力學的某齒輪軸沖擊應力分析[J].振動、測試與診斷,2013,33(2):273-276.

Yu Yinxin, Yuan Huiqun, Liang Mingxuan, et al. Impact stress analysis of gear shaft based on explicit dynamic[J]. Journal of Vibration, Measurement & Diagnosi, 2013,33(2):273-276. (in Chinese)

[8] 秦大同,周志剛,楊軍,等.隨機風載作用下風力發電機齒輪傳動系統動態可靠性分析[J].機械工程學報,2012,48(3):1-8.

Qin Datong, Zhou Zhigang, Yang Jun, et al. Time-dependent reliability analysis of gear transmission system of wind turbine under stochastic wind load[J]. Journal of Mechanical Engineering, 2012,48(3):1-8. (in Chinese)

[9] 張慶偉, 張博,王建宏,等.風力發電機齒輪傳動系統的動態優化設計[J].重慶大學學報,2010,33(3):30-35.

Zhang Qingwei, Zhang Bo, Wang Jianhong, et al. Dynamic optimization design of gear transmission system for wind turbine[J]. Journal of Chongqing University, 2010,33(3):30-35. (in Chinese)

[10]Maki K, Sbragio R, Vlahopoulos N. System design of a wind turbine using a multi-level optimization approach[J]. Renewable Energy, 2012,43:101-110.

[11]Chaari F, Baccar W, Abbes M S, et al. Effect of spalling or tooth breakage on gearmesh stiffness and dynamic response of a one-stage spur gear transmission[J]. European Journal of Mechanics, 2008,27(4):691-705.

[12]Sainsot P, Velex P, Duverger O. Contribution of gear body to tooth deflections-a new bidimensional analytical formula[J]. Journal of Mechanisms, 2004,126(4):748-752.

[13]Yang D C H, Sun Z S. A rotary model for spur gear dynamics[J]. Journal of Mechanisms Design, 1985,107(4):529-535.

[14]付才高.轉子動力學及整機振動[M]. 北京:航空工業出版社,2000:67-68.

10.16450/j.cnki.issn.1004-6801.2015.06.004

*國家自然科學基金資助項目(51165019)；甘肅省自然科學基金資助項目(1308RJYA018)

2013-11-16；

2014-03-03

TH132.413; TH113.1

劉宏,男，1986年4月生，博士研究生。主要研究方向為風力發電機機械系統動力學、多體系統動力學分析。 E-mail：hongliu2006@163．com