一種風電機組軸承簡化建模方法的研究

0 引言

滾動軸承的應用非常廣泛，不僅涉及到航天、船舶等領域，也涉及到風力發(fā)電設備[1]。對滾動軸承的分析起初采用的是動力學方法，對軸承的整個運動過程建立動力學方程，求得任意瞬間軸承各部件的運動參數(shù)和接觸應力。目前，軸承的分析主要采用計算機建模，通過仿真的形式分析軸承的各相關性能指標[2]，但滾動軸承的仿真建模主要還是以實體建模為主，所需計算量非常大。因此，本文針對風電機組上的變槳軸承與偏航軸承，采用GAP單元簡化軸承建模的方法，并用理論方法驗證其正確性。結果表明，這種簡化方法可有效模擬滾動軸承受力特性，對結構尺寸較大的軸承建?？善鸬揭欢ǖ慕梃b作用。

1 采用MSC.Marc/Mentat軟件中的GAP單元對軸承滾動體進行分析

四點接觸軸承有3個關鍵點，分別是滾動體中心點、滾動體與外圈接觸點、滾動體與內圈接觸點。四點接觸軸承在受力的過程中，其荷載傳遞形式是直線形式，如圖1所示。

圖1 荷載傳遞路線

根據(jù)荷載的傳遞路徑特征，在對軸承進行簡化建模時，需充分考慮MSC.Marc/Mentat中的GAP單元特性，通過對單元特性的設置，完成最終模型的建立。

GAP單元支持MSC.Marc/Mentat軟件中的任意其他單元連接，分為左右2個連接點，分別連接在模型的不同節(jié)點上，用以描述模型中任意兩節(jié)點的相對關系(包含摩擦、滑動或熱影響等)[3]。

1.1 實體滾動軸承參數(shù)及三維模型

采用三維軟件對滾動軸承(以某一單排四點接觸球軸承為例)進行建模，包含軸承內圈、軸承外圈、滾動體(模型中省略保持架，在后期的有限元分析軟件中采用一個微小彈簧代替保持架)來實現(xiàn)滾動體的固定。此模型建立比較簡單，具體如圖2所示。軸承各參數(shù)如表1所示。

圖2 四點接觸球軸承三維模型

表1 四點接觸球軸承參數(shù)表

1.2 滾動軸承的簡化建模分析

1.2.1 GAP單元的特性

假設B點代表外圈接觸點，C點代表內圈接觸點，那么B、C兩點可組成GAP單元，B、C兩點之間的距離可設定為滾動體的直徑Dw。當軸承受外界作用力時，GAP單元可出現(xiàn)荷載－位移關系，其關系可通過接觸面的切向分量和法向分量來表述[4]。在接觸區(qū)域，為了防止接觸體相互侵入，GAP單元的法向剛度將變得足夠大；而在未接觸區(qū)，它對分析對象的運動狀態(tài)不產生影響[5]。

當接觸點對的相對位移大于初始間隙時，法向力Fn的值就會是零，而接觸點對就會處于分離狀態(tài)，這時GAP單元不會發(fā)生力的傳遞，從而不會影響分析對象的運動狀態(tài)；而當相對位移小于初始間隙時，表明接觸點對處于接觸狀態(tài)，此時GAP單元就像一個線性彈簧，其法向接觸剛度為KA，在GAP單元的法線方向就會存在一個法向力Fn(負值)，內外圈就會通過GAP單元在B、C兩點之間傳遞荷載。

1.2.2 GAP單元對軸承的簡化建模

軸承三維模型建好之后，要對模型進行網格劃分。網格劃分是否合理將關系到后期模型計算的準確性和效率，所以網格要過渡合理、均勻；同時網格類型的選擇也至關重要，因網格類型有四面體、六面體、八面體等，每種類型適用于不同的模型，此處是對單個軸承進行分析，為了達到使用較少的網格重劃分次數(shù)而得到較高計算精度的目的，最終選擇六面體網格。

因軸承整圈具有左右對稱性，只取半個軸承模型進行分析，滾動體也只取半個進行實體建模(只考慮此滾動體的應力應變情況，與理論計算進行對比)，其他的滾動體只考慮力的傳遞，不再具體分析每個滾動體的應力應變，在這里將其他滾動體用“GAP單元”來替代定義。剛度是滾動軸承的一項重要指標，而GAP單元就像一個線性彈簧，將其長度設置為滾動體直徑后，在軟件中再將其初始軸向剛度和徑向剛度設置為零(初始狀態(tài)不受力)，一旦軸承承受外界作用力，GAP單元可出現(xiàn)荷載－位移關系，軸向剛度和徑向剛度就會出現(xiàn)一定的數(shù)值，接觸剛度這時就可近似模擬滾動體剛度[6]。

接觸角是滾動軸承的一個主要設計參數(shù)，在分析軸承受力、變形、運動關系和確定軸承承載能力時，都需要首先確定接觸角的大小。GAP單元由2個節(jié)點連接組成，一定數(shù)量的GAP單元可以近似代替一個滾動體，也就是在建立之初就已經指定了力的傳遞點，所以GAP單元無法精確描述接觸角的具體數(shù)值。這里針對接觸角不做具體分析，但GAP單元越多，越能更準確地反映滾動體和滾道之間的接觸情況。

GAP單元可以模擬任意直徑滾動體建模，只需在軟件中將GAP單元長度設置為滾動體直徑即可。

所謂一定數(shù)量的GAP單元代替一個滾動體，并不是說GAP單元越少越好，也不是越多越好。GAP單元太少時，不足以模擬滾動體與滾道的接觸狀態(tài)；而GAP單元太多時，就會增加軟件計算時間，降低效率。筆者前期已做過大量實驗驗證，為了使簡化模型更好地與理論吻合，GAP單元需要均勻布置到接觸角范圍內，相鄰GAP單元之間的角度差應控制在1°～3°之間。本模型采用14對GAP單元代替1個滾動體，相鄰GAP單元之間的角度差約為2.14°(本模型軸承接觸角為 30°)。

網格劃分單元的大小根據(jù)位置有不同的取值，其在內外圈接觸的地方大小約為2 mm，其他地方單元格可適當放大，如圖3所示。

圖3 基于GAP單元的滾動軸承簡化建模

1.2.3 保持架的模擬

在MSC.Marc/Mentat軟件中，選用剛度很小的彈簧單元對滾動體進行固定，小彈簧的作用類似于保持架，將各滾動體固定，并且均勻隔開。為了準確模擬保持架的作用，彈簧的一端與滾動體上節(jié)點相連，另一端設置為接地不動，如圖4所示。

圖4 彈簧單元固定滾動體

1.2.4 接觸設置

接觸設置主要是滾動體與內、外圈之間的狀態(tài)設定。首先，針對GAP單元，根據(jù)軟件中的要求，不需要對其進行接觸設置，只需將其與內、外圈節(jié)點連接即可；其次，對半個滾動體與內、外圈的接觸設置，這里設置為接觸，即“touching”狀態(tài)；最后，對其他參數(shù)進行設置，因為內圈、外圈和滾動體都假定采用軸承鋼，材料為GCr15，所以將它們設置為可變形體，相應的摩擦系數(shù)設置為0.15。

1.2.5 邊界條件、荷載的設置

1)邊界條件的設置：滾動軸承受力時，其運動狀態(tài)是不固定的，可以在軸承外圈上設置約束條件，同時對軸承內圈施加荷載；也可以在軸承內圈上設置約束條件，同時對軸承外圈施加荷載。在這里采取對外圈進行固定的方式，即對外圈采用全約束方式，同時約束軸承剖面沿法向的平動自由度，以及約束內圈沿軸向的轉動自由度。

2)荷載的設置：以上邊界條件設置好之后，開始在軸承內圈上加載來計算軸承接觸部位各應力變化。根據(jù)MSC.Marc/Mentat軟件對MPC點的介紹，這種點可以作為一個荷載施加中心，通過與模型相連，將荷載傳遞到模型上。此例中，將MPC點設置在軸承中心點處，并且將此點與內圈上表面所有點相連，荷載的傳遞路徑就是通過這些連線進行傳遞。MPC點設置好之后，在其上面施加混合荷載(包含翻轉力矩、徑向力和軸向力)。

邊界條件的設置及荷載的施加如圖5及表2所示。

圖5 邊界條件、荷載設置示意圖

表2 荷載值大小

1.2.6 滾動軸承簡化建模有限元計算結果

根據(jù)表2中的相關參數(shù)，可以得到滾動軸承的應力和應變結果，其中，內圈與滾動體之間的最大接觸應力值約為4603.2 MPa，最大位移值約為0.1511 mm；外圈與滾動體之間的最大接觸應力值約為4521.6 MPa，最大位移值約為0.1465 mm。

內、外圈與滾動體接觸的應力云圖如圖6、圖7所示。

圖6 內圈與滾動體接觸的應力云圖

圖7 外圈與滾動體接觸的應力云圖

2 滾動軸承理論計算

理論計算采用的軸承模型與有限元計算模型一致，參數(shù)表和荷載表也采取表1、表2中的數(shù)據(jù)。

2.1 徑向位移δr、軸向位移δa和傾斜角θ數(shù)值的確定

四點接觸軸承在受力時，會發(fā)生一定的位移。要通過理論方法計算軸承在荷載作用下的相關量，首先需要確定3個未知量：徑向位移δr、軸向位移δa和傾斜角θ，具體如圖8所示。

圖8 3個未知量示意圖

針對以上3個未知量，需要根據(jù)文獻[7-9]對其進行如下推導求解。圖9是在混合荷載作用下，對軸承外圈固定，軸承內圈曲率中心軌跡的變化曲線。其中，S為受載時內、外圈滾道曲率中心之間的最短距離。圖10為滾動體與滾道未加載與受載后的接觸示意圖，初始接觸角α0的取值為π/4。

圖9 內、外圈滾道曲率中心改變后軌跡示意圖

圖9、圖10中的其他主要參數(shù)可通過下文中的公式進行確定：

式中，δi為受載后內圈曲率中心位移量；δo為受載后外圈曲率中心位移量。

所以，受載后內、外圈曲率中心變化總量δn為：

圖10 滾動體與滾道未加載與受載后的接觸示意圖

內圈與滾動體之間的主曲率之和∑ρi可表示為：

軸承內圈輔助變量F()ρi可表示為：

圖11 內圈與滾動體接觸示意圖

式中，ρ1為滾動體與內圈接觸時的主曲率，ρ2為滾動體與外圈接觸時的主曲率，且為內圈滾道圓弧曲率，為軸承內圈溝底曲率，

根據(jù)式(3)、式(4)的推導，同理可得出軸承外圈輔助變量F()ρo及外圈與滾動體之間的主曲率之和

內圈滾道曲率中心點半徑Ri為：

外圈滾道曲率中心點半徑Ro：

滾動體總共有84個，假設其中某個滾動體的中心在X軸上，設其與X軸夾角為ψ，則初始值ψ0=0，其他滾動體在滾道圓周上均勻分布，那么繞其逆時針方向分布的第i個滾動體，其的范圍在[0,2π]之間，如圖 12 所示。

圖12 滾動體在滾道上的分布

根據(jù)參考文獻[9]，圖9中S值可表示為：

將式(7)帶入式(2)可得：

某個滾動體上所受荷載Qψi為：

式中，K為滾動體與滾道接觸的總剛度值；對于球軸承而言，n取1.5。

K可由式(10)求得：

式中，Ki為軸承內圈與滾動體之間的相對剛度；Ko為軸承外圈與滾動體之間的相對剛度。

其中：

對每個滾動體上的力進行分解，使每個滾動體上的受力均分解為與初始徑向力Fr0、軸向力Fa0以及力矩M0相同方向的分力，然后對各分力在各自方向上分別求和，具體如下：

式中，α為軸承受載后的接觸角，則有：

由此可推導出式(18)～式(20)：

軸承所受荷載如表2所示，根據(jù)式(16)～式(18)，可分別計算軸承的初始徑向位移δr′(僅受徑向力時)、初始軸向位移δa′(僅受軸向力時)、初始傾斜角θ′(僅受彎矩時)。

計算出各初始參數(shù)后，可根據(jù)牛頓迭代法進行收斂計算，具體如下[10]：

其中：

根據(jù)以上計算所得的3個初始值，將式(16)～式(18)聯(lián)立，采用牛頓迭代法進行求解。求解的過程是迭代收斂的過程，當最終收斂后，就會得到軸承在外荷載(見表2)作用下的徑向位移δr、軸向位移δa和傾斜角θ。

根據(jù)表1、表2中的數(shù)值及以上相關公式，可求得軸承其他參數(shù)值，具體如表3所示。

表3 計算所得參數(shù)值大小

2.2 接觸應力、位移求解

首先計算84個滾動體所受的荷載大小，可通過計算所得的軸承3個未知量(軸向位移δa、徑向位移δr和傾斜角θ)，結合公式(9)得到。針對以上有限元的計算結果，理論計算為了能與有限元模型更好的匹配，選擇同一位置處的滾動體進行荷載理論計算，即取ψ=π，則計算所得滾動體的荷載Q=58708.58383。

2.2.1 軸承內圈與滾動體接觸的最大應力及位移值

根據(jù)式(3)、式(4)及表 3可知，∑ρi=0.061145401，F(xiàn)(ρi)=0.6309204。

內圈與滾動體接觸的最大接觸應力pimax為：

根據(jù)F(ρi)的值，查《球軸承的設計計算》[9]中的赫茲接觸系數(shù)表，并利用插值法可得：πeaiebi= 1.96082836 ×10-3、eδi=2.48944154× 10-4，則pimax=4614.5 MPa。

滾動體與內圈的最大位移δi為：

將相關數(shù)值代入式(24)可得δi=0.1481 mm。

2.2.2 軸承外圈與滾動體接觸的最大應力及位移值

參照軸承內圈與滾動體的接觸計算，同理可計 算 出 ∑ρo=0.05948664、F(ρo)=0.6206288、則軸承外圈與滾動體接觸的最大接觸應力值最大位移δo=0.1475 mm。

3 對比分析

基于理論計算的結果與GAP單元對軸承的簡化建模計算結果進行對比，具體如表4所示。

表4 結果對比

根據(jù)表4中的數(shù)值可知，GAP單元建模計算的結果與理論計算的結果誤差較小，基本吻合。

4 結論

本文采用有限元分析方法，用GAP單元模擬軸承內、外圈的接觸，然后在添加外荷載的情況下對模型進行靜強度分析計算；同時，采用理論方法，在同種工況下對該軸承分析計算，針對兩種不同方法的計算結果分析對比，可得到以下結論：

1)GAP單元建模計算結果與理論計算結果誤差較小，吻合度高。

2)GAP單元模擬軸承滾動體進行建模，可實現(xiàn)軸承內、外圈力的傳遞。

3)當不針對滾動體進行分析時，可用GAP單元代替所有滾動體進行簡化建模，不會影響其他連接部件的計算結果。

4)GAP單元無法精確模擬軸承接觸角的具體數(shù)值，只能近似模擬滾動體與滾道間的接觸情況。要研究接觸角時，需采取實體建模進行分析。

5)GAP單元模擬滾動體對軸承的簡化建模，克服了實體軸承建模難度大、復雜性高的缺點，同時也大幅提高了計算速度。

6)這種簡化方法可以應用于風力發(fā)電機組中變槳軸承及偏航軸承的建模，從而為風電軸承建模提供一種新的思路。