張建宇 ，李加軍 ，杜曉鐘 ，胥永剛

（1.北京工業(yè)大學(xué)先進(jìn)制造技術(shù)北京市重點(diǎn)實(shí)驗(yàn)室，北京 100124；2.北京工業(yè)大學(xué)北京市精密測(cè)控技術(shù)與儀器工程技術(shù)研究中心，北京 100124；3.太原科技大學(xué)機(jī)械工程學(xué)院，山西 太原 030024）

1 引言

由于我國(guó)面臨的環(huán)境污染和能源枯竭問題日益嚴(yán)峻，作為一種清潔能源技術(shù)，風(fēng)力發(fā)電近年來獲得了長(zhǎng)足發(fā)展。其中，雙饋型機(jī)組由于發(fā)電效率高，因而贏得了更多的市場(chǎng)青睞。但是，隨著機(jī)組功率不斷提高，機(jī)組結(jié)構(gòu)愈發(fā)龐大，結(jié)構(gòu)內(nèi)部的耦合效應(yīng)更為突出，加之變速變載、極端風(fēng)況等因素的影響，機(jī)組內(nèi)部的振動(dòng)問題日益凸顯，由此誘發(fā)的傳動(dòng)鏈故障對(duì)機(jī)組安全運(yùn)行構(gòu)成了嚴(yán)重威脅。因此，如何準(zhǔn)確描述風(fēng)電機(jī)組內(nèi)部，特別是主傳動(dòng)系統(tǒng)的動(dòng)力響應(yīng)特征，是當(dāng)下的一個(gè)研究熱點(diǎn)。

長(zhǎng)期以來，眾多的國(guó)內(nèi)外學(xué)者都對(duì)風(fēng)力發(fā)電機(jī)的動(dòng)力學(xué)行為給予了重點(diǎn)關(guān)注。文獻(xiàn)[1]為了研究怎樣減小傳動(dòng)鏈各部件的機(jī)械應(yīng)力，建立了一個(gè)風(fēng)力機(jī)傳動(dòng)鏈的兩質(zhì)量塊模型，并采用仿真和試驗(yàn)結(jié)果驗(yàn)證了模型的正確性。文獻(xiàn)[2]建立了三種不同的風(fēng)力機(jī)行星齒輪箱的模型（包括剛性模型、柔性多體模型和純扭轉(zhuǎn)模型），利用有限元和實(shí)驗(yàn)方法對(duì)其進(jìn)行了研究。文獻(xiàn)[3]采用集中質(zhì)量法建立了風(fēng)力機(jī)傳動(dòng)鏈耦合動(dòng)力學(xué)模型，模型中考慮了柔性支撐，通過仿真計(jì)算獲得了風(fēng)力機(jī)各齒輪的振動(dòng)響應(yīng)特征。文獻(xiàn)[4]首先求解了風(fēng)力機(jī)齒輪傳動(dòng)系統(tǒng)在隨機(jī)風(fēng)作用下的動(dòng)力學(xué)特性，然后對(duì)其可靠性進(jìn)行了分析，但是建立的模型結(jié)構(gòu)較簡(jiǎn)單，只能對(duì)實(shí)際情況進(jìn)行近似求解和預(yù)測(cè)。文獻(xiàn)[5]建立風(fēng)力機(jī)傳動(dòng)鏈有限元模型，對(duì)其進(jìn)行動(dòng)力學(xué)特性研究，并對(duì)承載較大的斜齒輪部分進(jìn)行了柔性體分析，得到了關(guān)鍵部件的應(yīng)力云圖分布，為風(fēng)力機(jī)傳動(dòng)系統(tǒng)的疲勞壽命計(jì)算提供了參考。文獻(xiàn)[6]以單自由度系統(tǒng)簡(jiǎn)化齒輪傳動(dòng)，建立了8 自由度的風(fēng)力機(jī)傳動(dòng)鏈扭轉(zhuǎn)振動(dòng)模型，并完成了固有特性分析。文獻(xiàn)[7]考慮了齒輪副的時(shí)變嚙合剛度、傳動(dòng)誤差以及齒側(cè)間隙等內(nèi)部因素對(duì)系統(tǒng)的影響，分析了行星輪系的非線性動(dòng)力學(xué)特性。文獻(xiàn)[8]采用自回歸方法建立了隨機(jī)風(fēng)場(chǎng)的風(fēng)速模型，求得風(fēng)速時(shí)程曲線，以時(shí)變風(fēng)速引起的葉輪轉(zhuǎn)矩作為系統(tǒng)的輸入激勵(lì)，研究了變載荷激勵(lì)下的動(dòng)態(tài)響應(yīng)特性。

以某1.5MW 雙饋型風(fēng)力發(fā)電機(jī)為研究對(duì)象，建立整機(jī)的平移-扭轉(zhuǎn)耦合動(dòng)力學(xué)模型，重點(diǎn)考慮主傳動(dòng)鏈各機(jī)械部件的剛度激勵(lì)作用，以及時(shí)變風(fēng)載引起的外部激勵(lì)，研究各傳動(dòng)齒輪的扭轉(zhuǎn)振動(dòng)響應(yīng)特征，為機(jī)組的動(dòng)力特性評(píng)估提供基礎(chǔ)數(shù)據(jù)。

2 主傳動(dòng)系統(tǒng)的動(dòng)力學(xué)建模

研究的對(duì)象是1.5MW 風(fēng)電機(jī)組的主傳動(dòng)系統(tǒng)，其結(jié)構(gòu)包含一級(jí)行星輪系和兩級(jí)平行輪系，機(jī)組工作參數(shù)及傳動(dòng)結(jié)構(gòu)參數(shù)，如表1～表2 所示。該機(jī)組的結(jié)構(gòu)示意圖，如圖1 所示。圖中：b—葉輪；pi（i=1～3）—行星輪；r—內(nèi)齒圈；s—太陽(yáng)輪；1、2、3、4—兩級(jí)平行軸齒輪；g—發(fā)電機(jī)。

表1 風(fēng)電機(jī)組的工作參數(shù)Tab.1 Working Parameters of the Wind Turbine

表2 風(fēng)力機(jī)主傳動(dòng)系統(tǒng)幾何參數(shù)Tab.2 Geometric Parameters of the Main Drive System of the Wind Turbine

圖1 風(fēng)力機(jī)主傳動(dòng)系統(tǒng)示意圖Fig.1 Schematic Diagram of Main Drive System of Wind Turbine

2.1 建模的基本思路

針對(duì)的傳動(dòng)結(jié)構(gòu)，如圖1 所示。文獻(xiàn)[9]完成了建模與仿真分析，但該模型將齒輪2、3 等效為一個(gè)齒輪，而且沒有考慮嚙合阻尼的影響。為了全面描述系統(tǒng)的動(dòng)力學(xué)特性，文中建立了一個(gè)21自由度模型，示意圖，如圖2 所示。除考慮各齒輪對(duì)的嚙合剛度、軸系扭轉(zhuǎn)剛度、軸承支撐剛度以外，齒輪副的嚙合阻尼也計(jì)入其中。

圖2 風(fēng)力機(jī)主傳動(dòng)系統(tǒng)動(dòng)力學(xué)模型Fig.2 Dynamic Model of Main Drive System of Wind Turbine

模型（圖2）中的參數(shù)詳述如下：kj為軸承的支撐剛度（j=1，2…6）；kzj為軸的扭轉(zhuǎn)剛度（j=1，2，3，4）；ksi（t）、kri（t）、k12（t）、k34（t）分別為太陽(yáng)輪與行星輪、內(nèi)齒圈與行星輪、齒輪1 與齒輪2、齒輪3與齒輪 4 的時(shí)變嚙合剛度；csi、cri、c12、c34分別為太陽(yáng)輪與行星輪、內(nèi)齒圈與行星輪、齒輪1 與齒輪2、齒輪3 與齒輪4 的嚙合阻尼。

2.2 動(dòng)力學(xué)方程的建立

建立的動(dòng)力學(xué)模型包含9 個(gè)旋轉(zhuǎn)自由度和12 個(gè)平動(dòng)自由度，其中重點(diǎn)研究各級(jí)傳動(dòng)齒輪的振動(dòng)特性，因此，首先建立沿嚙合線方向的齒輪彈性變形表達(dá)式：

式中：δri，δsi，δ12，δ34—內(nèi)齒圈與行星輪、太陽(yáng)輪與行星輪、齒輪 1 與齒輪2、齒輪3 與齒輪4 沿嚙合作用線方向的彈性變形；

ηn—行星輪的位置角；α—齒輪壓力角。

根據(jù)達(dá)朗伯原理，作用在同一結(jié)構(gòu)上的慣性力、阻尼力、彈性力與外部激勵(lì)構(gòu)成平衡力系，因此，獲得各傳動(dòng)部件的動(dòng)力學(xué)方程如下：

（1）結(jié)構(gòu)件（包括葉輪和行星架）的振動(dòng)微分方程

（2）行星輪的振動(dòng)微分方程

其中，式（7）為扭轉(zhuǎn)振動(dòng)方程，式（8）和式（9）為垂直振動(dòng)方程。

（3）太陽(yáng)輪的振動(dòng)微分方程

其中，式（10）為太陽(yáng)輪扭轉(zhuǎn)振動(dòng)方程，式（11）和式（12）為垂直振動(dòng)方程。

（4）齒輪1 的振動(dòng)微分方程

其中，式（13）為齒輪 1 的扭轉(zhuǎn)振動(dòng)方程，式（14）和式（15）為垂直振動(dòng)方程。

（5）齒輪2 的振動(dòng)微分方程

其中，式（16）為齒輪 2 的扭轉(zhuǎn)振動(dòng)方程，式（17）和（18）為垂直振動(dòng)方程。

（6）齒輪3 的振動(dòng)微分方程

其中，式（19）為齒輪 3 的扭轉(zhuǎn)振動(dòng)方程，式（20）和（21）為垂直振動(dòng)方程。

（7）齒輪4 的振動(dòng)微分方程

其中，式（22）為齒輪 4 的扭轉(zhuǎn)振動(dòng)方程，式（23）和（24）為垂直振動(dòng)方程。

（8）發(fā)電機(jī)運(yùn)動(dòng)微分方程

式中：Jj—轉(zhuǎn)動(dòng)慣量（j=b，c，pi，s，1，2，3，4，g）；Φj—扭振角位移（j=b，c，pi，s，1，2，3，4，g）；mj—質(zhì)量（j=pi，s，1，2，3，4）；rj—基圓半徑（j=c，pi，s，1，2，3，4）；α—齒輪壓力角；zj、yj—在 z、y 方向上產(chǎn)生的振動(dòng)位移（j=pi，s，1，2，3，4）。

3 主傳動(dòng)系統(tǒng)的激勵(lì)載荷建模

風(fēng)力發(fā)電機(jī)的動(dòng)力來源于外部的時(shí)變風(fēng)載，通過葉輪轉(zhuǎn)化為主軸承受的時(shí)變扭矩。對(duì)風(fēng)力機(jī)的主傳動(dòng)系統(tǒng)來說，除了來自葉輪的扭矩激勵(lì)，系統(tǒng)內(nèi)部的剛度激勵(lì)，尤其是齒輪的時(shí)變嚙合剛度，也是結(jié)構(gòu)振動(dòng)的主要根源。因此，為了完成系統(tǒng)的響應(yīng)分析，必須建立準(zhǔn)確的激勵(lì)載荷模型。

3.1 內(nèi)部激勵(lì)建模

風(fēng)力機(jī)主傳動(dòng)系統(tǒng)的內(nèi)部激勵(lì)主要包括各齒對(duì)的時(shí)變嚙合剛度/嚙合阻尼、各軸承的支撐剛度和各軸系的扭轉(zhuǎn)剛度等。

3.1.1 齒輪副的時(shí)變嚙合剛度

目前，嚙合剛度的計(jì)算方法已經(jīng)比較成熟，常用方法包括能量法、石川法和有限單元法等。但是，風(fēng)力發(fā)電機(jī)包含一級(jí)行星輪系和兩級(jí)平行輪系，內(nèi)部多達(dá)8 個(gè)接觸對(duì)，如沿用上述方法逐一求解將過于繁瑣。因此，為了提高計(jì)算效率，節(jié)約計(jì)算內(nèi)存，根據(jù)嚙合剛度的周期特征，選用傅里葉級(jí)數(shù)求解各齒對(duì)的嚙合剛度，計(jì)算方法，如式（26）[9]所示。其中，kav—各齒輪副的平均嚙合剛度；ε—齒輪重合度；f—各齒對(duì)的嚙合頻率，數(shù)值大小，如表3 所示。顯然，式中的n 取值越大，計(jì)算值越接近實(shí)際的嚙合剛度（對(duì)直齒輪尤其如此），這里n 取值1000。

表3 齒輪副嚙合頻率Tab.3 The Meshing Frequency of Gear Pair

表中：fpr、fps、f12、f34—行星輪與內(nèi)齒圈、行星輪與太陽(yáng)輪、齒輪 1 與齒輪2、齒輪3 與齒輪4 的嚙合頻率。

3.1.2 齒輪副的嚙合阻尼

齒輪副的嚙合阻尼選用式（27）[8]進(jìn)行求解：

式中：ξ—嚙合阻尼比，一般?。?.03～0.17）；rp、rg—主、從動(dòng)齒輪的基圓半徑；Jp、Jg—主、從動(dòng)齒輪的轉(zhuǎn)動(dòng)慣量。

3.1.3 軸承的支撐剛度和軸系扭轉(zhuǎn)剛度

僅考慮軸承的支撐特性和軸的扭矩傳遞作用，因此，其剛度值均可視為常數(shù)，計(jì)算方法沿用剛度計(jì)算的經(jīng)驗(yàn)公式以及材料力學(xué)的經(jīng)典理論[11]結(jié)果，如表4 所示（參數(shù)含義見前述說明）。

表4 各軸及軸承的剛度計(jì)算結(jié)果Tab.4 Stiffness of Each Shaft and Rolling Bearings

3.2 外部激勵(lì)建模

自然界的風(fēng)瞬息萬變，風(fēng)速變化與地形、地貌、季節(jié)等多種因素有關(guān)。目前，模擬風(fēng)速的模型有很多，例如Weibull 分布、瑞利分布[12]概率統(tǒng)計(jì)模型。為了更好地說明風(fēng)速構(gòu)成的物理含義，選用組合風(fēng)速[9]求解風(fēng)力發(fā)電機(jī)承受的風(fēng)載荷，該模型將風(fēng)速變化視為平均風(fēng)、陣風(fēng)、漸變風(fēng)和隨機(jī)風(fēng)疊加的結(jié)果，其中既考慮了風(fēng)速中的直流分量，又考慮了風(fēng)速的隨機(jī)特性，因此，物理含義更為明確。

3.2.1 平均風(fēng)速

平均風(fēng)速是某一地區(qū)某一時(shí)段風(fēng)速的綜合反映，其數(shù)學(xué)模型可表示為：

3.2.2 陣風(fēng)風(fēng)速

陣風(fēng)是描述風(fēng)速突然變化的特性，其數(shù)學(xué)模型：

式中：vz—陣風(fēng)風(fēng)速；maxvz陣風(fēng)風(fēng)速的最大值；T1—陣風(fēng)的啟動(dòng)時(shí)間；T2—陣風(fēng)周期。

3.2.3 漸變風(fēng)風(fēng)速

漸變風(fēng)是描述風(fēng)速隨時(shí)間緩慢變化的特性，其數(shù)學(xué)模型：

式中：vj—漸變風(fēng)風(fēng)速；maxvj—漸變風(fēng)風(fēng)速的最大值；T′—漸變風(fēng)保持時(shí)間—漸變風(fēng)開始變化時(shí)間—漸變風(fēng)結(jié)束變化時(shí)間。

3.2.4 隨機(jī)風(fēng)風(fēng)速

隨機(jī)風(fēng)反映風(fēng)速變化的隨機(jī)性，用隨機(jī)噪聲的模型模擬，其數(shù)學(xué)模型可表示為：

式中：vs—隨機(jī)風(fēng)風(fēng)速；maxvs—隨機(jī)風(fēng)風(fēng)速的最大值；φv—（0～2π）之間均勻分布的隨機(jī)數(shù)；Ram（-1，1）—1 和-1 之間均勻分布的隨機(jī)數(shù)；ωv—風(fēng)速波動(dòng)的頻率，一般取（0.5π～2π）。

最后，將上述四種風(fēng)速成分進(jìn)行疊加，即可獲得組合風(fēng)速的模型。

根據(jù)某風(fēng)場(chǎng)風(fēng)速數(shù)據(jù)的統(tǒng)計(jì)結(jié)果，采用組合風(fēng)速數(shù)據(jù)，如表5 所示。僅研究風(fēng)力發(fā)電機(jī)在額定風(fēng)速下的振動(dòng)響應(yīng)特性，因此式（28）中平均風(fēng)速取值13m/s。根據(jù)組合風(fēng)速模型，完成60s 的風(fēng)速模擬，結(jié)果，如圖3 所示。

表5 組合風(fēng)速計(jì)算參數(shù)Tab.5 Calculation Parameters of Combined Wind Speed

圖3 時(shí)變風(fēng)速曲線Fig.3 The Curve of Time-Varying Wind Speed

根據(jù)空氣動(dòng)力學(xué)的相關(guān)理論，時(shí)變風(fēng)速與葉輪扭矩之間存在一定的換算關(guān)系，其計(jì)算模型，如式（33）所示?；谠撚?jì)算式和上述風(fēng)速曲線，可得葉輪傳遞給風(fēng)力機(jī)主軸的外部扭矩，如圖4所示。

式中：Tin—葉輪扭矩；Ct（λ）—扭矩系數(shù)；v—風(fēng)速；R—風(fēng)力機(jī)葉輪半徑。

圖4 主傳動(dòng)系統(tǒng)外部扭矩曲線Fig.4 External Torque of Main Drive Systemin Time Domain

4 系統(tǒng)的仿真結(jié)果與分析

基于前述動(dòng)力學(xué)模型，以及內(nèi)、外部激勵(lì)載荷的分析結(jié)果，采用四階-五階Runge-Kutta 方法仿真研究了風(fēng)力機(jī)主傳動(dòng)系沿各方向的動(dòng)態(tài)響應(yīng)特性。囿于篇幅限制，此處僅討論各齒輪的扭振響應(yīng)。

在外部風(fēng)載和內(nèi)部激勵(lì)的聯(lián)合作用下，行星輪系構(gòu)件（即行星輪和太陽(yáng)輪）的扭轉(zhuǎn)振動(dòng)角位移，如圖5 所示?？梢?，兩者的整體振動(dòng)趨勢(shì)受外部激勵(lì)載荷（即圖4 所示扭矩曲線）的影響；而響應(yīng)波形中包含的周期特征，則來源于輪系的內(nèi)部激勵(lì)，換言之，圖中的振動(dòng)變化周期與行星輪-太陽(yáng)輪的嚙合剛度時(shí)變周期一致。

圖5 行星輪系振動(dòng)響應(yīng)波形圖Fig.5 Vibration Response of Planetary Gear Train in Time Domain

兩級(jí)平行輪系中各齒輪組件的扭振響應(yīng)位移，其中，蘊(yùn)含的規(guī)律與行星輪系相同，即整體趨勢(shì)反映外部激勵(lì)的作用，周期特征與內(nèi)部激勵(lì)的周期相對(duì)應(yīng)，如圖6 所示。

圖6 平行輪系的振動(dòng)響應(yīng)波形圖Fig.6 Vibration Response of Parallel Wheel System in Time Domain

此外，如果將圖5 與圖6 的響應(yīng)幅值進(jìn)行對(duì)比，可以發(fā)現(xiàn)：行星輪系的扭振響應(yīng)明顯小于平行輪系。換句話說，系統(tǒng)振動(dòng)由低速級(jí)向高速級(jí)逐漸增大，其中以高速軸齒輪（即齒輪4）的幅值為最大。為了進(jìn)一步揭示各級(jí)齒輪的振動(dòng)特性，提取了齒輪構(gòu)件的角加速度響應(yīng)結(jié)果，并計(jì)算其響應(yīng)譜結(jié)果，如圖7～圖12 所示。

圖7 行星輪角加速度振動(dòng)響應(yīng)Fig.7 Angular Acceleration Vibration Response of Planetary Gear

圖8 太陽(yáng)輪角加速度振動(dòng)響應(yīng)Fig.8 Angular Acceleration Vibration Response of Solar Gear

圖9 齒輪1 角加速度振動(dòng)響應(yīng)Fig.9 Angular Acceleration Vibration Response of Gear 1

圖10 齒輪2 角加速度振動(dòng)響應(yīng)Fig.10 Angular Acceleration Vibration Response of Gear 2

圖11 齒輪3 角加速度振動(dòng)響應(yīng)Fig.11 Angular Acceleration Vibration Responseof Gear 3

圖12 齒輪4 角加速度振動(dòng)響應(yīng)Fig.12 Angular Acceleration Vibration Response of Gear 4

由各輪系的扭轉(zhuǎn)加速度響應(yīng)結(jié)果，可以發(fā)現(xiàn)：位移響應(yīng)中的振動(dòng)趨勢(shì)線（與外部激勵(lì)相關(guān)），經(jīng)過兩次求導(dǎo)運(yùn)算，在加速度響應(yīng)中得以消除；另一方面，主要受外部激勵(lì)影響的振動(dòng)幅值，仍然呈現(xiàn)由低速級(jí)向高速級(jí)逐漸增大的趨勢(shì)，統(tǒng)計(jì)各齒輪角加速度響應(yīng)的有效值結(jié)果，如表6 所示。振動(dòng)能量的變化趨勢(shì)更為明確。

在各組角加速度波形圖中，周期振動(dòng)特征都異常明顯，進(jìn)一步分析響應(yīng)譜，如圖7（b）～圖12（b）所示?？梢园l(fā)現(xiàn) 25.2Hz、83.8Hz和326.5Hz 等比較突出的頻率成分，這些數(shù)值分別與行星輪-太陽(yáng)輪、齒輪1-齒輪2、齒輪3-齒輪4 的嚙合頻率相一致，體現(xiàn)了各對(duì)齒輪嚙合剛度的激勵(lì)作用，表明仿真結(jié)果是真實(shí)、可信的。

表6 動(dòng)角加速度有效值統(tǒng)計(jì)（rad/s2）Tab.6 RMS Values of Angular Acceleration Response

5 結(jié)論

通過建立1.5MW 風(fēng)力發(fā)電機(jī)的平移-扭轉(zhuǎn)耦合動(dòng)力學(xué)模型，同步施加外部風(fēng)載和內(nèi)部剛度/阻尼的聯(lián)合激勵(lì)作用，仿真獲得各輪系的扭轉(zhuǎn)振動(dòng)響應(yīng)特征，結(jié)論如下：

（1）在各齒輪構(gòu)件的扭振位移響應(yīng)中，振動(dòng)趨勢(shì)及響應(yīng)幅值主要受外部激勵(lì)風(fēng)載的影響，而內(nèi)部的齒輪嚙合剛度激勵(lì)則與振動(dòng)響應(yīng)的周期特征有關(guān)。

（2）在扭轉(zhuǎn)角加速度響應(yīng)中，振動(dòng)趨勢(shì)線得以消除，但振動(dòng)幅值仍能體現(xiàn)外部激勵(lì)的影響；加速度響應(yīng)的振動(dòng)周期與施加的剛度激勵(lì)周期嚴(yán)格對(duì)應(yīng)，齒輪的嚙合頻率在各組響應(yīng)譜中均較為突出。

（3）對(duì)比行星輪系與兩級(jí)平行輪系的仿真結(jié)果，可以發(fā)現(xiàn)：無論是位移響應(yīng)，還是角加速度響應(yīng)，均呈現(xiàn)由低速級(jí)向高速級(jí)逐漸增大的趨勢(shì)，其中高速軸齒輪的響應(yīng)幅值最突出，據(jù)此可判斷該齒輪發(fā)生局部破壞的可能性較高。

（4）僅提供了風(fēng)力機(jī)主傳動(dòng)系的扭振響應(yīng)結(jié)果，今后將結(jié)合系統(tǒng)的垂振響應(yīng)分析，全面評(píng)估結(jié)構(gòu)的振動(dòng)特征和劣化可能性。