基于無網格法的風力發電機主軸強度仿真研究*

聶新宇,汪小芳,劉 勇,白聰兒,孫哲杰,夏國俊

(1.浙江運達風電股份有限公司,浙江 杭州 310012;2.浙江省風力發電技術重點實驗室,浙江 杭州 310012;3.浙江大學 能源工程學院,浙江 杭州 310027)

0 引 言

隨著碳達峰、碳中和戰略目標的提出,風電行業又將迎來快速發展的良好機遇[1]。為了降低度電成本,風力發電機組的單機容量被不斷提高,風輪直徑與載荷也相應增大[2-3]。

作為機組傳動鏈的關鍵一環,主軸負責將源自風輪的扭矩傳遞至齒輪箱、發電機,將傾覆力矩、推力等其他載荷傳遞至主軸承等支撐結構[4-5]。這些載荷的增大會顯著提高主軸不同區域的應力,因此,為了保障機組的安全運行,基于強度仿真的主軸結構設計顯得非常重要。

目前,關于風力發電機主軸的強度分析主要以有限元法為主。

王斌等人[4]利用ANSYS建立了主軸疲勞強度的分析模型,探究了卸荷槽寬度及深度等尺寸對鎖緊螺紋段疲勞損傷的影響。汪亞洲等人[6]建立了主軸的極限強度有限元分析模型,并借此對軸承前沿大圓角和橢圓長徑尺寸進行了優化,降低了主軸的最大應力。趙佰余等人[7]采用Workbench計算了主軸的極限與疲勞強度,并對可設計區進行了結構優化,在主軸減重的同時還提高了其安全余量。白儒等人[8]通過有限元分析,從結構剛度協調性和表面粗糙度兩方面對主軸進行了優化,在滿足強度要求的前提下,實現了主軸的輕量化設計目標。

為準確模擬載荷的傳遞路徑,這些主軸模型中均包含主軸承結構,通常為雙列調心滾子軸承(SRB)或單列圓錐滾子軸承(TRB),包含數十至上百個滾子,一般采用桿單元、彈簧單元或間隙單元對其作等效處理[9-10],設定僅受壓屬性。其建模與計算過程較為復雜耗時。

為了提高仿真效率,近年來,基于無網格法的仿真方式在電子設備、車輛、醫學等領域得到越來越多的應用[11]。例如,MAGALHES F C等人[12-13]對3D打印、器官等復雜結構進行了靜強度評估,簡化了網格劃分過程,有效識別了結構集中應力。顧張祺等人[14]對電子產品組件開展了模態分析與加速度過載分析,簡化了幾何清理、接觸設置等環節,有效縮短了產品研制初期的設計時間。

由上述研究可知,無網格仿真可以簡化傳統結構仿真中的前處理環節,適合對含接觸的大型裝配體強度進行快速評估。然而,目前該方法在含軸承的主軸結構上還少有應用。

基于上述原因,筆者將以單SRB和雙TRB兩類機型的主軸結構為例,分別建立其有限元和無網格仿真模型,探究適用于三維原型的無網格模型前處理方法,分析無網格法在主軸裝配體上的計算可靠性,及其在建模過程與計算時間上的優勢,進而提出一種利用無網格法提高迭代效率的主軸結構設計方法。

1 主軸有限元模型

1.1 模型概述

筆者利用ANSYS軟件建立主軸的有限元模型,采用GL2010規定的輪轂轉動坐標系[15]。

輪轂轉動坐標系如圖1所示。

圖1 輪轂轉動坐標系

圖1中,坐標系隨輪轂轉動,原點在輪轂中心,XR為轉動軸方向,ZR朝向葉片,并與XR垂直,YR垂直于XR和ZR,根據右手定則確定其方向。

兩類機型的主軸裝配結構如圖2所示。

圖2 兩類主軸裝配結構示意圖

裝配體一般包含主軸、輪轂、軸承、擋圈、脹緊套及齒輪箱接口法蘭等部件:單SRB機型只含一個主軸承,配合齒輪箱扭力臂支撐主軸;雙TRB機型主要通過2個軸承對主軸進行支撐[10]。

ANSYS網格劃分圖如圖3所示。

圖3 ANSYS網格劃分圖

主軸結構:為準確模擬其集中應力,筆者采用高階六面體單元SOLID186劃分網格;輪轂結構:做了適當簡化并采用較疏的網格。

SRB結構模型單元總數為368 343,節點總數為1 481 364;TRB結構單元總數為431 014,節點總數為1 638 017。

1.2 材料屬性

主軸是整體鍛件,材料為合金結構鋼42CrMo4,彈性模量為2.1×105MPa,泊松比為0.3;

輪轂材料為QT400,彈性模量為1.69×105MPa,泊松比為0.275;

其他部件采用鋼結構,材料屬性參照主軸設定。

1.3 軸承等效

軸承滾子數量較多,且與內外滾道間為非線性接觸。為了減少計算量,筆者采用多組受壓不受拉的桿單元LINK10等效模擬[7],通過調整桿單元的截面積,使其剛度和滾子、滾道間的接觸剛度一致。

根據ISO 16281[16],滾子與滾道間為線接觸,其彈性常數CL可通過下式計算得到:

(1)

式中:Lwe為滾子有效接觸長度。

進而可得接觸載荷Q與位移δ間的關系式為:

Q=CLδ10/9

(2)

等效后的軸承模型如圖4所示。

圖4 軸承有限元模型剖面圖

1.4 邊界條件

模型邊界條件設置如下:擋圈環向面與主軸、軸承間采用標準接觸,摩擦系數為0.2,其他相鄰部件間均為綁定接觸;軸承安裝于軸承座內,外側環向面可視為固定約束。在實際運行中,主軸后端與齒輪箱相連,其軸向轉動受到限制。

因此,筆者在齒輪箱重心與脹緊套或齒輪箱接口法蘭間建立剛性連接,約束重心點ROTX=UY=UZ=0;風輪受到重力和風動載荷等,利用Bladed軟件計算其各向力矩達到極限時的載荷,包括Mx_max、My_max、Mz_max和Myz_max4種工況;在輪轂中心建立加載點,通過剛性梁單元與變槳軸承的安裝面相連,將載荷施加于該加載點。

兩類主軸的4種極限工況載荷如表1所示。

表1 極限載荷表

2 無網格模型

2.1 背景介紹

在有限元仿真中,存在網格劃分復雜、大型裝配體計算耗時的問題。近年來,一些簡化甚至略去網格劃分的數值分析方法被提了出來,工程上統稱為無網格法。

無網格法的理論基礎多樣:如外部有限元逼近法(external finite element approximation,EFEA),采用式(3)所示近似函數,在單元邊界自由度基礎上附加內部自由度,可用不規則幾何體離散求解域,而非四面體、六面體等標準單元,因此較容易采用軟件實現自動網格劃分[17];此外還有隱式邊界法(implicit boundary method,IBM),采用不受分析對象結構限制的規則網格,對完全處于結構內部、部分處于結構內部的網格,采取不同算法進行分析計算[18]。

筆者將以基于外部有限元逼近法的SimSolid軟件為例,進行主軸的建模分析。

外部有限元逼近法的近似函數如下:

(3)

式中:uh為近似解;ai為單元內部自由度;pi為內部自由度基函數;bk為單元邊界自由度;pk為邊界自由度基函數。

2.2 建模方法

SimSolid建模步驟如下:1)導入幾何體;2)設定材料參數;3)創建接觸關系;4)建立分析模型。

材料、載荷等參數與上文有限元模型一致;所不同的是:1)不含網格劃分過程;2)導入模型為幾何體原型,主軸承不做簡化。

在該無網格模型中,接觸是通過幾何體間的接觸對實現。接觸對的建立由如下3個參數控制:間隙容差、穿透容差和連接分辨率。合理的參數設置要保證滾子與滾道間充分接觸,同時剛度不至太高。

模型邊界條件設置方法如下:

1)遠程載荷施加。以輪轂中心為遠程點,輪轂上的3個變槳軸承的安裝面為加載面,施加表1所示載荷。

2)主軸尾端約束。在脹緊套或齒輪箱接口法蘭端面上,建立“接地襯套”形式的虛擬連接,設定軸向轉動剛度與垂直軸向的平動剛度為較大數值(此處設定為1020N·m/rad和1020N/m)。

3)軸承移動約束。在軸承外側環向面上施加鉸鏈約束,保證軸承可以繞軸向轉動,同時軸向、徑向平動自由度為零。

以SRB模型為例,滾子與滾道間的接觸對呈長條狀。軸承內部接觸示意圖如圖5所示。

圖5 軸承內部接觸示意圖

3 結果與討論

3.1 剛度分析

有限元分析是業內公認的主軸強度分析方法[4,7,10]。因此,筆者以有限元結果為基準,判斷無網格模型的準確度。

首先,比較結構變形。此處以Myz_max工況為例,兩類結構在不同仿真方式下的計算合位移,如圖6所示。

由圖6可以看出:對于單SRB結構,SimSolid計算位移分布規律與ANSYS一致,其最大值均出現在輪轂端部,分別是31.06 mm和31.61 mm,相對誤差為1.77%,該處遠離軸承和齒輪箱的約束,因此位移較大;

對于雙TRB結構,最大位移同樣出現在輪轂端部,數值分別是12.06 mm和11.88 mm,相對誤差為-1.49%。

在不同工況下,2種模型的最大合位移對比如表2所示。

表2 最大合位移對比表

由表2可知:最大合位移的相對誤差在-6.06%～2.28%之間。在仿真模型中,輪轂與主軸綁定,主軸的變形受到軸承的影響,軸承剛度越大,對主軸的變形限制越明顯,輪轂位移量越小;反之增大。

位移相近,說明無網格法對軸承剛度的模擬較為可靠。

3.2 單SRB結構集中應力分析

在剛度計算準確的基礎上,筆者分析了無網格法計算主軸應力的準確度。

此處,筆者首先以單SRB結構Myz_max工況為例。

2種模型計算得到的主軸Mises應力分布如圖7所示。

圖7 單SRB結構Myz_max工況下主軸Mises應力云圖

由圖7可以看出:主軸應力均呈對稱分布,對稱軸為合彎矩Myz方向,主軸首尾兩端應力較高。前端圓角位于軸承與輪轂之間,受彎矩最大,因此應力水平普遍較高。疊加該區域的截面積變化梯度大,應力集中現象較為明顯。主軸尾端半徑最小,受扭矩引起的切應力最高。

筆者統計這兩個區域不同工況下的最大Mises應力結果,如圖8所示。

圖8 單SRB結構主軸不同工況下2種模型應力對比圖

由圖8可以看出:2種模型的最大Mises應力相對誤差在-2.36%～2.94%之間,最大值均出現在Myz_max工況下的尾端,分別是295 MPa和289 MPa;

主軸所用合金鋼的屈服強度σs=560 MPa,根據GL2010規定,材料安全系數γm取1.1[15],由式(4)可得主軸許用應力[σ]=509 MPa,2種模型計算主軸應力均滿足極限強度要求。

許用應力計算表達式如下:

[σ]=σs/γm

(4)

3.3 雙TRB結構集中應力分析

筆者接下來分析雙TRB結構主軸的應力。

同樣以Myz_max工況為例,ANSYS與SimSolid模型計算得到的主軸Mises應力分布,如圖9所示。

由圖9可以清楚地看到:主軸最大應力均位于前端圓角處,分別是147 MPa和148 MPa,相對誤差為0.68%,尾端最大應力分別是90 MPa與92 MPa,相對誤差為2.22%(如表1所示);該工況下,扭矩在各向載荷中的占比較低,導致尾端應力要明顯小于前端圓角。

筆者統計了該結構的2種模型在不同工況下的應力計算結果,如圖10所示。

圖10 雙TRB結構主軸不同工況下2種模型應力對比圖

由圖10可以看出:2種模型最大Mises應力相對誤差在-6.71%～4.76%之間,最大值均出現在My_max工況下的前端圓角處,分別是168 MPa和176 MPa,同樣滿足極限強度設計要求。

3.4 計算效率對比

在建模步驟、計算時間與結果方面,筆者對2種模型的差異進行了統計,結果如表3所示。

表3 2種模型對比表

由表3可知:相較于ANSYS,基于無網格法的SimSolid建模省去了網格劃分和軸承等效等步驟,大大縮短了建模時間,降低了建模難度;

在計算時間上,針對單SRB和雙TRB結構,ANSYS計算時間平均為13 min 44 s和29 min 24 s,SimSolid為3 min 12 s和3 min 36 s,相對用時縮短了77%與88%,有限元法計算時間受模型規模、結構復雜度影響較大,無網格法則相對較小,展現了其在風力發電機組結構越來越復雜背景下的應用潛力。

3.5 無網格輔助主軸結構設計方法

鑒于無網格法在計算效率上的優勢,筆者提出一種無網格輔助主軸結構設計方法。

該主軸結構設計方法的流程圖如圖11所示。

圖11 無網格輔助主軸結構設計流程圖

在新型尤其是大功率主軸的研發過程中,基于經驗設計的初始結構往往和最終結構差異較大,需經歷多次迭代過程,以得到滿足極限、疲勞強度要求的結構。

由于有限元分析的迭代成本較高,因此,筆者在傳統的設計流程中添加圖11虛線框所示無網格法,利用無網格法計算高效的特點,預先優化主軸結構,減少后續迭代次數,從而縮短整體設計時間[19]。

4 結束語

無網格仿真可以簡化傳統結構仿真中的前處理環節,然而,目前該方法在含軸承的主軸結構上還少有應用。

為了提高風力發電機中傳遞載荷的關鍵部件—主軸的設計效率,筆者以單SRB和雙TRB兩類機型的主軸結構為例,分別建立了有限元仿真模型和無網格仿真模型,分析無網格法在主軸裝配體上的計算可靠性,及其在建模過程與計算時間上的優勢,進而提出了一種利用無網格法提高迭代效率的主軸結構設計方法。

研究結論如下:

1)對于兩類機型,無網格法計算主軸變形結果與有限元分析結果相吻合,最大位移相對誤差在-6.06%～2.28%之間,說明無網格法在主軸及相關軸承的剛度模擬上的準確度較高;

2)比較2種模型可知,在主軸前端圓角和尾端等應力集中區域,兩類機型的最大應力相對誤差分別在-2.36%～2.94%及-6.71%～4.76%之間,說明無網格法應力計算結果是可靠的;

3)在保證計算精度的同時,無網格法的建模過程更加簡便,同時計算時間顯著減低,由此,提出了一種無網格輔助主軸結構設計的方法。

在后續的研究中,筆者擬將無網格法應用于輪轂、機架等含復雜曲面,以及裝配關系更為復雜的結構件中,以擴展該方法的應用范圍,進一步提升風力發電機組的結構設計效率。