汲柏良，杜昱成，秦清海

(曲阜師范大學工學院，山東 日照 276800)

風電作為綠色清潔能源，是全球實現(xiàn)低碳、可持續(xù)發(fā)展的有效途徑之一。我國的低風速區(qū)域分布廣，占總面積的近7成[1]。近年來，國內低風速區(qū)的風電裝機量迅速增加，進一步開發(fā)低風速區(qū)風能資源已經成為風電行業(yè)的共識。

風電機組按照風力機與發(fā)電機的傳動鏈劃分，可分為直驅、半直驅以及雙饋3種類型[2]。半直驅方案采用“齒輪箱-永磁發(fā)電機-變頻器”的技術路線，吸納了直驅型與雙饋型的長處，既在一定程度上避免了雙饋型使用多級增速齒輪箱導致的故障率高的問題，又通過提高發(fā)電機輸入轉速相對降低了發(fā)電機的體積、重量及制造成本，半直驅漸漸成為兼顧生產制造成本與發(fā)電效益的最優(yōu)解。

齒輪箱是一種工作于高空、低速環(huán)境下，受不規(guī)律風載提高發(fā)電機轉速的機械傳動裝置。齒輪箱故障會造成風電機組停機，降低供電可靠性，對風電行星齒輪傳動系統(tǒng)進行動力學分析具有重要的現(xiàn)實意義。Kahraman最早求解了包含齒側間隙、齒形誤差和時變嚙合剛度等因素的數(shù)學模型表達形式[3]，建立了行星齒輪純扭轉動力學模型并求解分析[4]。文獻[5-6]中推導了行星齒輪的平移-扭轉動力學模型，文獻[6]考慮了齒根裂紋對行星齒輪的影響，并對其進行了分析。文獻[7]以人字齒行星齒輪為研究對象，推導了動力學微分方程模型，將受嚙合相位影響的時變嚙合剛度考慮其中，并分析受不均勻載荷影響時的變化規(guī)律。齒輪的嚙合運動中，嚙合的齒對數(shù)及嚙合的位置隨時間發(fā)生變化，因而齒輪的嚙合剛度具備時變特征[8]。

本文主要研究低風速半直驅風電行星齒輪箱。首先對滿足10 kW半直驅風力發(fā)電機運行條件下單級行星齒輪的各主要構件進行參數(shù)匹配與結構設計，而后基于經典力學定律，計入時變嚙合剛度與綜合嚙合誤差的影響，通過集中參數(shù)模型推導行星齒輪的平移-扭轉時變非線性動力學微分方程組，最后通過四階Runge-Kutta法對動力學微分方程組進行求解分析。

1 行星齒輪結構設計

半直驅風電機組如圖1所示，主要包含三部分：第一部分為風力機，是風能捕獲裝置；第二部分為風電齒輪箱，是增速裝置；第三部分為永磁同步發(fā)電機，是機電能量轉換裝置。

圖1 半直驅風電機組

該風電機組主要運行于低風速工況下，輸入轉速低，轉矩大。2Z-X型直行星齒輪結構簡單、承載力強、效率高，因此將2Z-X型直行星齒輪作為風電齒輪箱的增速機構。行星齒輪所需的性能要求如表1所示。

表1 風電行星齒輪性能要求

(1)

式中：Kd為算式系數(shù)，取值為768；T1為嚙合齒輪副中小齒輪的名義轉矩；KA為使用系數(shù)，與隨機載荷有關；KH∑為綜合系數(shù)；KHP為載荷分布系數(shù)，與接觸強度有關；φd為小齒輪齒寬系數(shù)；σHlim為接觸強度下齒輪的接觸疲勞極限；u為齒數(shù)比，“+”表示外嚙合，“-”表示內嚙合。

(2)

式中：Km為算式系數(shù)，取值為12.1；KFΣ為綜合系數(shù)；KFP為載荷分布系數(shù)，與彎曲強度有關；YFa1為小齒輪齒形系數(shù)；z1為齒輪副中小齒輪齒數(shù)；σFlim為彎曲強度下齒輪的接觸疲勞極限。

選取適當?shù)膮?shù)分別代入式(1)與式(2)中計算可得：m1=4.33；m2=4.55。根據國標(GB/T 1357—2008)選定模數(shù)為5 mm。經計算，風電行星齒輪的系統(tǒng)參數(shù)如表2所示。

表2 風電行星齒輪系統(tǒng)參數(shù)

2 行星齒輪動力學模型

2.1 動力學建模

該風電行星齒輪傳動系統(tǒng)中，內齒圈保持靜止，扭矩由行星架側輸入軸輸入，通過帶動行星輪由太陽輪側輸出軸輸出。

為簡化分析，在建立平移-扭轉耦合模型時忽略輸入軸與輸出軸對行星齒輪的影響；另行星齒輪為直齒輪，沒有軸向力的參與，不計行星齒輪的軸向振動；3個行星齒輪的參數(shù)完全一致，并且是均載的；不計齒側間隙以及摩擦對行星齒輪動態(tài)響應的作用；動力學模型中的阻尼都假定為線性阻尼。

基于上述假設，建立行星齒輪的平移-扭轉耦合動力學模型如圖2所示。圖中有3個坐標系：坐標系OXY是靜止坐標系；Oxy是旋轉坐標系，其原點固定，兩坐標軸與行星架同向等速旋轉；onxnyn是以行星輪中心點為坐標原點，兩坐標軸與坐標系Oxy兩坐標軸分別平行且與行星輪等角速度旋轉的動坐標系，n=1，2，3。

圖2 平移-扭轉耦合動力學模型

行星齒輪各構件包含3個自由度，2個平動x、y與1個扭轉u。太陽輪、行星輪以及齒圈的直徑分別取值為各自基圓的直徑，行星架基圓的半徑取值為太陽輪中心點與行星輪中心點間的距離。行星齒輪各構件經彈簧和阻尼器的并聯(lián)結構相連接。ksp、csp分別表示太陽輪與行星輪間的嚙合剛度與嚙合阻尼；krp、crp分別表示齒圈與行星輪間的嚙合剛度與嚙合阻尼(p=1、2、3)。kl、cl(l=s、c、r、p)分別表示各個構件的支承剛度與支承阻尼；kcp、ccp表示行星架與行星輪間的支承剛度與支承阻尼。

以Δ表示彈性構件間的壓縮形變量，太陽輪與行星輪間及齒圈與行星輪間的壓縮形變量為

(3)

(4)

式中：αs、αr為太陽輪與行星輪及齒圈與行星輪間的嚙合角(按照標準中心距安裝時，嚙合角等于分度圓壓力角)；φi為各行星輪安裝相位角；es(t)、er(t)為太陽輪與行星輪及齒圈與行星輪間的綜合嚙合誤差。

第i個行星輪與太陽輪及第i個行星輪與齒圈間的嚙合力為

(5)

式中：ksi為第i個行星輪與太陽輪間的嚙合剛度；csi為第i個行星輪與太陽輪間的嚙合阻尼；kri為第i個行星輪與齒圈間的嚙合剛度；cri為第i個行星輪與齒圈間的嚙合阻尼。

第i個行星輪與行星架間的形變量為

(6)

式中：αi為第i個行星輪的壓力角。

行星輪與行星架間在x、y向上的支承力為

(7)

式中：kcpi為行星輪與行星架間的支承剛度；ccpi為行星輪與行星架間的支承阻尼。

各部件l在x、y及u向上的支承力為

(8)

式中：kl為部件l的支承剛度；cl為部件l的支承阻尼。

由此可得行星齒輪傳動系統(tǒng)動力學微分方程組見式(9)—(12)。

(9)

(10)

(11)

(12)

式中：ms、mc、mr、mp、Js、Jc、Jr、Jp、Ts、Tc、Tr、Tp、Rbs、Rbc、Rbr、Rbp分別為太陽輪、行星架、齒圈與行星輪的質量、轉動慣量、轉矩以及基圓半徑；ωc為行星架的轉速。

2.2 參數(shù)分析

2.2.1 時變嚙合剛度

太陽輪與行星輪間及齒圈與行星輪間，其嚙合剛度隨時間發(fā)生周期性變化，其主要原因如下。

a.嚙合時齒輪的嚙合對數(shù)發(fā)生改變。

b.輪齒嚙合過程中，其嚙合位置隨時間發(fā)生改變。

行星齒輪各構件間的嚙合剛度可使用與嚙合頻率有關的Fourier級數(shù)來等效，為簡化計算，此處取其前三次諧波。

(13)

式中：ksiav、kriav分別為行星輪與太陽輪及行星輪與齒圈間的平均嚙合剛度，可查閱GB/T 3480—1997計算得到；ksij、krij(j=1、2、3)為傅里葉系數(shù)；tsi、tri為太陽輪與各行星輪及齒圈與各行星輪嚙合的時間差異；fm為嚙合頻率，可由式(14)計算得到[11]。

(14)

式中：fc為行星架轉頻；Zr為齒圈齒數(shù)；Zs為太陽輪齒數(shù)；fs為太陽輪轉頻。

2.2.2 綜合嚙合誤差

綜合嚙合誤差產生于齒輪生產制造的過程中，在齒輪嚙合的過程中會造成沖擊響應。在動力學模型中，可采用正弦函數(shù)的形式來近似表示。

(15)

式中：Es為齒輪副綜合嚙合誤差的幅值；βs、βr為太陽輪與行星輪間及齒圈與行星輪間綜合嚙合誤差的初相位。

2.2.3 各構件質量及轉動慣量

行星齒輪中各構件質量可通過式(16)—(18)計算得到[14]。

(16)

式中：

dm=(da+df)/2

(17)

b=φdd1

(18)

式中：ρ為材料密度；d0為回轉軸直徑；B為齒寬；da為齒頂圓直徑；df為齒根圓直徑；φd為齒寬系數(shù)，d1為主動輪的分度圓直徑。

各個構件的轉動慣量為

(19)

3 計算結果分析

針對所選型設計的行星齒輪進行動力學分析，經計算各構件參數(shù)如表3所示。

表3 行星齒輪傳動系統(tǒng)參數(shù)

基于前文所推導的行星齒輪動力學模型，代入以上參數(shù)降階后采用四階五級定步長Runge-Kutta法編程求解。

3.1 確定工況下的動態(tài)響應分析

在不考慮負載的情況下，設定行星齒輪的輸入轉速為20 r/min，輸入轉矩為4775 N·m。行星齒輪各構件的動態(tài)響應分別如圖3(a)-(f)所示。

(a) 太陽輪x向振動位移

(b) 太陽輪y向振動位移

(c) 行星架x向振動位移

(d) 行星架y向振動位移

(e) 行星輪x向振動位移

(f) 行星輪y向振動位移圖3 行星齒輪各構件動態(tài)響應

由圖3(a)-(f)可知：各構件x向與y向的振動為非簡諧周期振動。各構件x向振動位移的時程曲線關于y=0近似對稱，太陽輪的x向振動位移幅值約為3.8 μm，行星架的x向振動位移幅值為0.96×10-3μm，行星輪的x向振動位移幅值為1.2×10-3μm。各構件y向振動位移中，太陽輪的y向振動位移幅值約為6.4 μm，行星架的y向振動位移幅值約為2.0×10-3μm，行星輪的y向振動位移幅值約為6.33 μm。

風電行星齒輪傳動系統(tǒng)往往運行在低速重載的環(huán)境下，系統(tǒng)中各個構件的振動幅度均在“微米”的數(shù)量級，可滿足要求。

3.2 輸入轉速的影響分析

風電行星齒輪的太陽輪通過輸出軸與發(fā)電機相連，太陽輪的動態(tài)響應關系到發(fā)電機的工況，故在不考慮負載條件下，設定輸入轉矩為4775 N·m，分別計算輸入轉速為10 r/min，20 r/min及30 r/min時的動態(tài)響應，得到太陽輪扭轉方向的振動位移時域如圖4所示，去除太陽輪扭轉方向振動位移的直流分量并分別進行頻域分析如圖5(a)-(c)所示。

圖4 不同轉速下太陽輪扭轉位移時域

(a) n=10 r/min

(b) n=20 r/min

(c) n=30 r/min圖5 不同輸入轉速下太陽輪扭轉位移頻譜

由圖4及圖5(a)-(c)可知：在保持輸入轉矩為定值4775 N·m時，改變行星齒輪的輸入轉速，其扭轉振動位移的峰值均為34.15 μm，由此可見定轉矩下，輸入轉速的變化對太陽輪扭轉位移幅值的影響非常??；行星齒輪中太陽輪扭轉位移的頻率分布與行星齒輪的嚙合頻率有關，主要為嚙合頻率及其倍頻，并且與轉速成線性關系。

3.3 輸入轉矩的影響分析

不考慮負載影響，設定輸入轉速為20 r/min，計算輸入轉矩分別為500 N·m、2000 N·m以及4775 N·m時太陽輪扭轉位移的動態(tài)響應，得到太陽輪扭轉方向的振動位移時域如圖6所示，同樣去除太陽輪扭轉方向振動位移的直流分量并分別進行頻域分析如圖7(a)-(c)所示。

圖6 不同轉矩下太陽輪扭轉位移時域

(a) T=500 N·m

(b) T=2000 N·m

(c) T=4775 N·m圖7 不同輸入轉矩下太陽輪扭轉位移頻譜

由圖6及圖7(a)-(c)可知：在輸入轉速不變情況下，當行星齒輪的輸入轉矩分別為500 N·m、2000 N·m以及4775 N·m時，太陽輪的扭轉位移幅值分別為3.57 μm、14.30 μm、34.15 μm，峰值分別為0.28 μm、1.13 μm、2.71 μm，因此太陽輪的扭轉位移會隨著輸入轉矩的增大而緩慢增大，扭轉位移的數(shù)量級仍處于“微米”級別；輸入轉矩僅對太陽輪扭轉位移的幅值產生正相關影響，與太陽輪扭轉位移頻率分布無關。

4 結語

a.風電行星齒輪在低速重載的工況下，太陽輪相比于齒圈與行星輪，其x向與y向振動位移相對較大，但其振動幅值仍維持在“微米”級。

b.在僅行星齒輪的輸入轉速發(fā)生變化的情況下，其扭轉位移的振動幅值不會發(fā)生變化，扭轉位移的頻率分布主要為嚙合頻率的倍頻，并且與轉速成強線性關系。

c.在行星齒輪的輸入轉矩發(fā)生變化的情況下，其扭轉位移的振動幅值隨著輸入轉矩而增大，扭轉位移的頻率分布并未隨轉矩的變化而發(fā)生改變。