風電機組齒輪箱傳動系統(tǒng)建模與動力學特性仿真分析

董 健，胡 亮，程志學，柳亦兵

（華北電力大學 能源動力與機械工程學院，北京 102206）

風電機組齒輪箱傳動系統(tǒng)建模與動力學特性仿真分析

董 健，胡 亮，程志學，柳亦兵

（華北電力大學 能源動力與機械工程學院，北京 102206）

針對目前風電機組常見的兩級行星輪系加一級平行軸傳動的齒輪箱結構，建立其傳動系統(tǒng)虛擬樣機模型?；贖ertz彈性接觸理論建立仿真齒輪嚙合剛度計算公式；借助ADAMS動力學分析軟件，對輪齒動態(tài)嚙合力進行仿真，研究輪齒嚙合力在齒輪轉動過程中的變化規(guī)律；將輸出級小齒輪柔性化，建立傳動系統(tǒng)的剛柔耦合模型，分析仿真過程中輸出級小齒輪的應力變化。利用ADAMS對復雜齒輪箱傳動系統(tǒng)進行動力學建模與仿真取得很好的效果，為傳動系統(tǒng)的動態(tài)設計及優(yōu)化提供依據。

振動與波；齒輪箱；虛擬樣機；動力學分析；動態(tài)響應

風電機組齒輪箱作為風電機組中的重要組成部分，其主要功能是傳遞扭矩和升速。行星齒輪傳動由于其結構緊湊、傳動比大、傳動平穩(wěn)等優(yōu)點，被廣泛用于機械傳動的各個領域。風電機組齒輪箱傳動系統(tǒng)多采用行星齒輪加平行軸齒輪復合傳動的形式。國內外學者對齒輪箱做了大量研究。李潤方、王建軍等在文獻[1]中系統(tǒng)的闡述了齒輪系統(tǒng)動力學的基本理論及分析方法。陳裴等[2]基于ADAMS多體動力學仿真軟件，對行星齒輪進行剛柔耦合仿真，分析了太陽輪發(fā)生故障后的動力學響應。朱才朝等[3]在考慮齒輪系統(tǒng)時變剛度、齒側間隙和制造誤差的基礎上，建立了具有多級齒輪傳動的大型風電齒輪箱的齒輪—傳動軸—箱體系統(tǒng)耦合非線性動力學模型，在考慮系統(tǒng)內部激勵的情況下對整個耦合系統(tǒng)動態(tài)特性進行了研究。Gang Liu，Robert G. Parker等[4]考慮齒輪副摩擦力及其彎曲效應，對單齒輪副的動態(tài)特性進行了分析。Woohyung Kim等在文獻[5]中分析了由軸承變形引起的齒輪時變壓力角及重合度對行星齒輪動態(tài)響應的影響。劉振皓、巫世晶等[6]以ADAMS軟件為平臺建立了復合行星齒輪傳動系統(tǒng)的虛擬樣機模型，對其動力學特性進行仿真研究，為改善系統(tǒng)動力學特性及實現系統(tǒng)的動態(tài)設計提供依據。李興然等[7]利用Pro/E軟件建立了配氣機構的三維實體模型，結合ANSYS軟件和ADAMS軟件，建立了配氣機構的多柔體動力學模型。但是目前對于兩級行星齒輪傳動系統(tǒng)動力學特性的研究文獻較少。

本文利用Pro/E和ADAMS軟件建立了兩級行星齒輪傳動加一級平行軸齒輪傳動系統(tǒng)的虛擬樣機模型，并對其動力學特性進行了仿真分析。

1 基于ADAMS的動力學理論

ADAMS中定義剛體的廣義坐標時采用質心的笛卡爾坐標和反應剛體方位的歐拉角，即考慮約束方程，ADAMS利用帶拉格朗日乘子的拉格朗日第一類方程的能量形式得到如下方程

T為系統(tǒng)廣義坐標表達的動能，qj為廣義坐標，Qj為廣義坐標qj方向的廣義力，最后一項為約束方程和拉格朗日乘子在廣義坐標qj方向的約束反力。

引入廣義動量Pj

由系統(tǒng)運動方程、約束方程構成的系統(tǒng)動力學微分代數方程為

其中Φ為約束函數，Φq為約束方程的Jacobian矩陣，λ為拉格朗日乘子，H為外力坐標轉換矩陣，F為系統(tǒng)所受外力。

為避免微分方程的Jacobian矩陣出現病態(tài)，可采用降階積分方法。利用SI 2降階積分方法對微分代數方程進行變換，得到公式

ADAMS中采用變系數的BDF剛性積分法對式（4）所示的動力學微分方程進行求解[8]，BDF剛性積分方法是一種自動變階、變步長的預估校正法，主要分為預估階段、校正階段和誤差控制階段。

在預估階段，根據當前時刻的系統(tǒng)狀態(tài)矢量值，用泰勒級數預估下一時刻系統(tǒng)的狀態(tài)矢量值

式中時間步長h=tn+1-tn。這種預估算法得到的新時刻的系統(tǒng)狀態(tài)矢量值通常不準確，可以由Gear k+ 1階積分求解程序來校正。

式中 yn+1為 y(t)在t=tn+1時的近似值；β0和 αi為Gear積分程序的系數值。式（6）經過整理，可以表示為

在校正階段，將狀態(tài)變量的值代入系統(tǒng)方程G(y,y·,t)=0，如果方程成立則狀態(tài)變量就為方程的解；若G(y,y·,t)≠0，利用修正的Newton-Raphson求解線性方程得到Δy，更新y，使系統(tǒng)方程更接近于成立。JΔy=G(y,y·,tn+1)，其中J為系統(tǒng)的雅克比矩陣。重復上述步驟，直到Δy足夠小。

在誤差控制階段，對狀態(tài)向量的預測值與校正值進行比較，如果誤差過大，則舍棄此解；如果誤差在允許范圍以內，則接受此解。不論誤差大小，都要對步長h和階數n進行求解，直到達到仿真結束時間停止。

2 傳動系統(tǒng)虛擬樣機模型

2.1 模型建立

本文所研究風電機組齒輪箱傳動系統(tǒng)為兩級行星齒輪加一級平行軸齒輪傳動形式，其結構簡圖如圖1所示，嚙合副1、3為內嚙合，嚙合副2、4、5為外嚙合。

傳動系統(tǒng)關鍵齒輪部件參數如表1所示，按表格參數利用三維繪圖軟件Pro/E建立傳動系統(tǒng)各部件三維模型并進行裝配，經檢驗確保系統(tǒng)無干涉后，另存為x_t文件格式導入到ADAMS中，定義材料屬性，得到傳動系統(tǒng)虛擬樣機模型如圖2所示。

2.2 添加約束

ADAMS約束庫中提供了多種約束類型，常用的約束有旋轉副（joint）、固定副（fixed）以及接觸（contact）。根據實際的運行情況在零件之間以及零件與大地（ground）之間添加約束：主軸與大地、太陽輪與大地、行星輪與行星架、輸出級齒輪與大地之間添加轉動副；行星輪與太陽輪、行星輪與齒圈、輸出級齒輪之間添加接觸；主軸和低速級行星架、低速級太陽輪和中間級行星架、齒圈和大地之間添加固定副。傳動系統(tǒng)虛擬樣機模型中共包含10個旋轉副、6個固定副、13個接觸副。

表1 傳動系統(tǒng)關鍵齒輪部件參數

圖1 傳動系統(tǒng)結構簡圖

圖2 傳動系統(tǒng)虛擬樣機模型

2.3 定義接觸

兩個構件的表面發(fā)生接觸時，構件就會在接觸的位置產生接觸力。ADAMS中可以通過定義接觸模擬輪齒的嚙合過程。ADAMS中計算接觸力的方法有沖擊函數法和補償法[9]。沖擊函數法根據Hertz彈性接觸理論，采用Impact函數提供的非線性等效彈簧阻尼模型作為接觸力的計算模型，接觸力由兩部分組成：一是由于兩構件相互滲入而產生的彈性力；另一個是由于相對速度產生的阻尼力。本文采用沖擊函數法計算輪齒接觸力?；跊_擊函數法的輪齒接觸力的計算公式為：

式中K為剛度系數；g為碰撞過程中相互滲入的深度；e為碰撞指數；C為阻尼系數；step為階躍函數；dmax為碰撞過程最大允許穿透深度。

剛度系數K取決于相互碰撞物體的材料和結構形式，計算公式為

式中R為綜合曲率半徑，E為綜合彈性模量，可由下式計算

將式（10）、式（11）代入到公式（9）得到

式中E1、E2為兩個接觸齒輪材料的彈性模量；μ1、μ2為兩個接觸齒輪材料的泊松比；R1、R2為兩接觸齒輪接觸點的當量半徑，由于齒輪的齒高和分度圓半徑相比較小，因此變動范圍不大，可近似以分度圓上的值來代替，這樣的近似誤差不大[10]。

齒輪的彈性模量E1=E2=2.1×105MPa，泊松比μ1=μ2=0.29。根據經驗，最大穿透深度dmax取0.2 mm，阻尼系數C取30N/(s·mm)，碰撞指數e取1.5。

2.4 斜齒輪受力分析

斜齒輪嚙合過程中，輪齒所受總法向力Fn處于與輪齒相垂直的法面上，可分解為軸向力Fa、徑向力Fr、圓周力Ft，計算公式為

式中αn為齒輪法向壓力角，β為齒輪螺旋角，T為齒輪所受轉矩，P為功率單位是kW，n為轉速，單位是r/min。

對行星齒輪與太陽輪嚙合副，其受力如圖3所示。

圖3 行星齒輪與太陽輪嚙合副受力分析

圖中Ts為太陽輪所受轉矩，ωs為太陽輪角速度，ω為行星架角速度，斜齒輪所受總法向力Fn在斜齒輪旋轉平面內的投影太陽輪所受輪齒接觸力沿y，z方向的接觸力為

3 動力學仿真

3.1 設置仿真條件

齒輪箱輸入轉速為17.4 r/min，總傳動比為104.125，結合各級齒輪齒數算出第一級行星齒輪嚙合頻率為29.58Hz，第二級行星齒輪嚙合頻率155.295 Hz，輸出級齒輪嚙合頻率為815.298 8 Hz。

為使輸入轉速不出現突變，在主軸上以STEP函數添加一個速度驅動在0.2 s時轉速達到104.4o/s即17.4 r/min，并保持不變。在輸出級齒輪添加一個負載轉矩大小為7.9×106N·mm，設置仿真時間為1 s，仿真步長為0.000 1 s。

3.2 運動學分析

圖4為齒輪箱傳動系統(tǒng)穩(wěn)定運行時各級齒輪角速度曲線，從仿真結果可以看出，穩(wěn)定運行時各級齒輪角速度在均值附近波動，包含有明顯的周期成分。

圖4 齒輪箱傳動系統(tǒng)各級齒輪角速度曲線

各級齒輪角速度仿真結果與理論計算結果如表2所示，仿真結果與理論計算結果十分接近，誤差很小。

表2 各級齒輪角速度仿真值與理論值對比（o/s）

3.3 動力學分析

圖5為低速級行星齒輪與太陽輪嚙合齒輪副作用在太陽輪上嚙合力矩的時域波形及其頻譜。從圖中可以看出傳動系統(tǒng)穩(wěn)定運行時，作用在太陽輪上的嚙合力矩在均值附近波動，并含有一定的周期成分，嚙合力矩的均值為5.21×107N·mm，表明漸開線齒輪在嚙合過程中存在沖擊振動。從頻譜中可以看出，嚙合力矩時域波形中主要頻率成分為第一級嚙合頻率的2倍頻，同時包含了第一級嚙合頻率1倍頻、3倍頻等成分。

圖5 低速級行星齒輪與太陽輪嚙合力矩時域波形及其頻譜

圖6為中間級行星齒輪與太陽輪嚙合齒輪副作用在太陽輪上的嚙合力矩的時域波形及其頻譜。從圖中可以看出，傳動系統(tǒng)穩(wěn)定運行時中間級太陽輪所受嚙合力矩同低速級太陽輪所受嚙合力矩具有同樣的規(guī)律：在均值上下波動，同時包含一定的頻率成分，嚙合力矩的均值為9.95×106N·mm。從其頻譜圖中可以看出，太陽輪所受嚙合力矩中包含有多種頻率成分，圖中譜線59.23 Hz對應第一級嚙合頻率的2倍頻，310.29 Hz對應第二級嚙合頻率的2倍頻。

輸出級嚙合齒輪副作用在輸出軸上的嚙合力矩的時域波形及其頻譜如圖7所示，從圖中可以看出作用在輸出軸上的嚙合力矩具有明顯的周期成分，嚙合力矩在一個均值上下波動，嚙合力矩的均值為8.024×106N·mm。從其頻譜中可以看出輸出軸所受嚙合力矩中明顯的頻率成分為30.54 Hz，對應輸出軸的轉頻，同時還包含了輸出級嚙合頻率的1倍頻以及2倍頻等成分。

圖6 中間級行星齒輪與太陽輪嚙合力矩時域波形及其頻譜

圖7 輸出級齒輪副嚙合力矩時域波形及其頻譜

表3為各級齒輪所受力矩的理論值與仿真值對比。從表中可以看出仿真得到各級齒輪所受嚙合力矩的均值與理論計算值十分接近。

表3 各級齒輪所受力矩仿真值與理論值對比（N·mm）

圖8為中間級太陽輪所受輪齒接觸力在y方向和z方向的分量。從圖中可看出中間級太陽輪所受y方向接觸力和z方向接觸力變化規(guī)律相似，相位差大約為90°，與理論分析結果一致。

圖8 中間級太陽輪所受接觸力在y、z方向分量

輪齒嚙合力中徑向力和圓周力分量對軸承的選型以及軸的優(yōu)化設計十分重要。圖9給出了輸出級齒輪副嚙合力圓周力分量的時域波形及其頻譜，圓周力分量的均值為-9.61×104N，從其頻譜圖中可以看出圓周力分量中主要的頻率成分為輸出級嚙合頻率的1倍頻成分，同時還包含低速級嚙合頻率的2倍頻、中間級嚙合頻率的2倍頻以及輸出級嚙合頻率的2倍頻成分等。

圖9 輸出級齒輪副嚙合力圓周力分量時域波形及其頻譜

輸出級齒輪副嚙合力徑向力分量的時域波形及其頻譜如圖10所示，徑向力分量均值為3.69×104N。對比徑向力分量頻譜圖和圓周力方向頻譜圖可以發(fā)現，信號中包含有相同的頻率成分，軸承選型及輸出軸的設計時應使軸承和軸的固有頻率避開圓周力和徑向力中包含的頻率成分。

圖10 輸出級齒輪副嚙合力徑向力分量時域波形及其頻譜

3.4 剛柔耦合分析

風電機組齒輪箱實際運行時，輸出級小齒輪發(fā)生故障的概率較高，因此將輸出級小齒輪柔性化后導入到ADAMS中建立傳動系統(tǒng)的剛柔耦合模型。模型仿真條件同剛體模型一致，仿真得到小齒輪最大應力為55.96 MPa。圖11為應力最大處節(jié)點的應力時間歷程曲線，從圖中可以看出應力變化具有一定的周期性。

圖11 應力最大處節(jié)點應力時間歷程曲線

4 結語

本文以風電機組齒輪箱中常見的兩級行星齒輪加一級平行軸齒輪傳動結構為研究對象，利用ADAMS動力學仿真軟件建立了其傳動系統(tǒng)的虛擬樣機模型，對傳動系統(tǒng)的動力學特性進行仿真分析，得出如下結論：

(1)仿真得到各級齒輪的角速度以及嚙合力幅值與理論計算值十分接近，誤差很??；

(2)穩(wěn)定運行時嚙合力中包含明顯的動載成分，從嚙合力頻譜中可以發(fā)現信號中包含了豐富的頻率成分，其中各級嚙合頻率的1倍頻、2倍頻等高頻成分突出，對傳動系統(tǒng)的的振動信號進行分析時應重點關注這些頻率成分；

(3)輸出級小齒輪危險點處的應力變化具有明顯的周期性，其最大應力小于小齒輪材料的屈服應力。

對常見的風電機組齒輪箱傳動系統(tǒng)虛擬樣機的動力學特性進行仿真研究，為改善傳動系統(tǒng)動力學特性，實現傳動系統(tǒng)的動態(tài)設計提供了一定的指導依據。

致謝

本文的工作得到中國華能集團公司科技項目“風機設備狀態(tài)監(jiān)測及故障診斷系統(tǒng)的研制”（HNKJ13-H20-05）資助，在此表示衷心的感謝。

[1]李潤方，王建軍.齒輪系統(tǒng)動力學[M].北京：科學出版社，1997.

[2]陳裴，荊建平.行星齒輪斷齒故障動力學仿真與故障提取[J].噪聲與振動控制，2014，34(05)：182-186.

[3]朱才朝，黃澤好，唐倩，等.風力發(fā)電齒輪箱系統(tǒng)耦合非線性動態(tài)特性的研究[J].機械工程學報，2005，41(8)：203-207.

[4]Gang Liu,Robert G Parker.Impact of tooth friction and its bending effect on gear dynamics[J].Journal of Sound and Vibration,2009,320(4-5):1039-1063.

[5]Woohyung Kim,Ji Yeong Lee,Jintai Chung.Dynamic analysis for a planetary gear with time-varying pressure angles and contact ratios[J].Journal of Sound and Vibration,2012,331(4):883-901.

[6]劉振皓，巫世晶，潛波，等.復合行星齒輪傳動系統(tǒng)虛擬樣機仿真研究[J].中國機械工程，2012，23(16)：1962-1966.

[7]李興然，張保成，續(xù)彥芳，等.配氣機構多柔體系統(tǒng)動力學仿真[J].噪聲與振動控制，2010，30(3)：74-77.

[8]陳立平.機械系統(tǒng)動力學分析及ADAMS應用教程[M].北京：清華大學出版社，2005.

[9]李增剛.ADAMS入門詳解與實例[M].北京：國防工業(yè)出版社，2012.

[10]龍凱，程穎.齒輪嚙合力仿真計算的參數選取研究[J].計算機仿真，2002，19(6)：16-22.

Modeling and Dynamic SimulationAnalysis of Wind Turbine Gearbox Transmission System

DONG Jian,HU Liang,CHENG Zhi-xue,LIU Yi-bing
（School of Energy,Power and Mechanical Engineering,North China Electric Power University, Beijing 102206,China)

Virtual prototype model of the drive system of a wind turbine generator with a common gearbox of twostage planetary and one-stage parallel shaft gear transmission was built.Based on the Hertz elastic contact theory,the formula of gear mesh stiffness was established.By means of the dynamics analysis software ADAMS,the dynamic performance of the gear meshing was simulated and the meshing force was obtained.Then,the rigid-flexible model of the transmission system with a flexible output-stage pinion was established,and the stress variation of the pinion was analyzed. It is shown that the method of using ADAMS for dynamic simulation and modeling of complicated gear transmission system is very effective.The simulation results can provide some useful reference for gearbox optimal design.

vibration and wave;gearbox;virtual prototype;dynamic analysis;dynamic response

TH113；TH132.4

A

10.3969/j.issn.1006-1335.2015.03.020

1006-1355(2015)03-0092-06

2014－12－12

董健（1982－），男，博士生，主要研究方向：機械系統(tǒng)動力學特性分析。E-mail:dongjian@gdupc.cn