隨著風力發(fā)電行業(yè)的迅速發(fā)展，葉片長度不斷增加，結構設計安全余量逐步降低。葉片作為核心的捕風裝置，其結構安全穩(wěn)定是風機正常運行的前提。根據葉片的結構特點，屈曲已成為葉片極限工況下的主要的失效形式[1,2]。由于全尺寸試驗檢驗葉片的承載能力與結構可靠性成本過高，故基于有限元的仿真分析一直是驗證葉片結構穩(wěn)定性的主要手段。

 

風力葉片研發(fā)伊始，研究者基于有限元方法針對葉片屈曲失效進行了大量的工作。薛彩虹等[3]基于 Abaqus 采用線性特征值屈曲方法分析了 40m 級葉片的穩(wěn)定性；李成友等[4]應用 Nastran軟件中特有的特征值抽取算法判斷 34m 級葉片的臨界失穩(wěn)點，并計算出屈曲因子；張立[5]等基于 ANSYS 軟件研究了 NREL 5MW 風力機葉片承受彎扭耦合載荷作用下的線性屈曲特征。不難發(fā)現(xiàn)，大部分研究都基于特征值屈曲的分析方法，特征值屈曲以小位移、小應變的線彈性理論為基礎，忽略了結構在受載變形過程中結構構形的變化[6,7]，而葉片結構復雜、鋪層繁多、變形很大，基于線性方法難免存在誤差，故 IEC61400-5 和 DNV-0376 均提出了基于非線性的葉片屈曲分析方法，以準確地預測葉片的承載能力，進而指導葉片結構鋪層設計[8]。在大兆瓦、小節(jié)圓、長葉片設計成為當前研發(fā)主流的情況下，非線性屈曲非常具有工程應用價值。

 

基于上述討論，本文應用 Abaqus 有限元軟件對時代新材 90m 級葉片進行屈曲分析,采用線性、非線性多種方法對大揮工況的臨界載荷進行求解，并對兩種方法比較分析。

 

2. 有限元模型

時代新材某 90m 級葉片為雙腹板結構，如圖 2-1 所示的葉片截面，蒙皮、主梁是葉片最主要的部分，主梁承擔葉片絕大部分的彎曲載荷由，而蒙皮起氣動作用且承擔部分彎曲載荷及大部分的剪切載荷，其中彎曲載荷大約為總載荷的 20% ～ 30%。

 

 

圖 2-1 90m 級葉片截面圖

 

2.1 葉片有限元模型

基于外部氣動外形的輸入，應用 Foucs 進行模型建立，依據鋪層設計劃分結構區(qū)域，并劃分網格，網格尺寸為 50mm*50mm，如圖 2-2 所示。將模型導入 Abaqus 有限元軟件中，應用python 與 abaqus 的接口進行編程并對葉片模型進行鋪層，其中應用的為 S4r 殼體單元以及C3D8R實體單元。

 

 

圖 2-2 90m 級葉片有限元模型

 

2.2 邊界條件

在實際運行中葉片與輪轂通過螺栓相連，結構分析需將其簡化為懸臂梁模型，故對葉根部進行固定全約束。對于葉片的載荷形式，根據整機氣動模擬給到的彎矩載荷等效轉換為集中力，施加在葉片主梁位置，如圖 2-3 所示。

 

 

圖 2-3 葉片有限元模型邊界條件

 

3. 屈曲分析方法

屈曲分析主要用于研究結構在特定載荷下的穩(wěn)定性，以及確定結構失穩(wěn)的臨界載荷，它包括線性屈曲和非線性屈曲分析。屈曲穩(wěn)定性分析的意義就是為進一步優(yōu)化結構提供依據。

 

3.1 特征值屈曲分析法

其中線性特征值屈曲分析法，它適用于對理想彈性結構的理想屈曲強度進行預測，主要是使用特征值公式計算造成結構負剛度的應力剛度陣的比例因子。結構在達到屈曲載荷之前其位移—變形曲線表現(xiàn)出線性關系，達到屈曲以后曲線將跟隨另外的路徑。發(fā)生轉折的這一點稱為分支點，分支點的載荷稱為屈曲載荷，結構發(fā)生屈曲的不同形態(tài)成為屈曲模態(tài)。

 

 

 

入即上述方程的最小特征值解，代表臨界載荷比例因子，最小特征值乘以初始載荷得到臨界載荷，即抗屈曲值。

 

3.1 非線性 Riks 分析法

Riks 是求解非線性屈曲應用較為廣泛的方法，其基本思想是假定載荷沿弧長方向增加或減小，Abaqus 中提供了修正的 Riks 算法，通過沿切線移動到給定距離到當前解點限制增量大小，進而結合獲得的點且與同一切線正交的平面中進行搜索取得平衡[7]，具體如下：

 

以Pⁿ表示加載模式，N表示模型的自由度，入載荷放大因子，故模型實時運算的載荷為

入Pⁿ為對應的位移，正則化為uⁿ，實時速度uo，同理，以矢量 (uⁿ；入)表達當前弧長步計算結果

 

 

圖 3-1 修正的 Riks 法位移-載荷示意曲線

 

 

 

其中， Abaqus/Standard 自動負載增量算法可根據收斂率針對靜態(tài)問題進行調整，之后，同一弧長步內，取得當前弧長步最后的解并進行迭代，進而得到最后結果。

 

4. 屈曲仿真結果與分析

以大揮方向的工況作為研究對象，首先基于 Abaqus 的 Buckle 求解整機廠輸入的載荷進行特征值屈曲分析，求解結果如圖 4-1 和圖 4-2所示，在這里列舉了一階模態(tài)及二階模態(tài)的屈曲線性計算結果，一階模態(tài)屈曲特征值為 1.978，二階模態(tài)屈曲特征值為 1.997，均在 1.965以上，滿足 GL2015 及 IEC61400-5 標準，也就證明葉片鋪層結構對該載荷結構相對穩(wěn)定，但值得注意的線性分析結果表面相對容易出現(xiàn)屈曲的位置靠近葉尖。

 

 

圖 4-1 葉片一階模態(tài)特征值屈曲

 

 

圖 4-2 葉片二階模態(tài)特征值屈曲

 

之后，應用非線線性弧長法進行計算，考慮到葉片由于生產制造過程中可能存在的誤差，引入幾何缺陷，應用 abaqus關鍵字語句*IMPERFECTION,FILE=Job-name,STEP=1 將特征值計算得到前兩階模態(tài)的節(jié)點擾動引入到修正的 Riks后屈曲計算中，依據葉片設計 IEC61400-5標準對缺陷進行縮放。非線性屈曲求解結果如下，當載荷放大因子為時間載荷的1.653倍時停止計算，圖 4-3 和 4-4 為計算結束后的應力及應變云圖。

 

 

圖 4-3 修正的 Riks 求解屈曲應力云圖

 

 

圖 4-4 修正的 Riks 求解屈曲應變云圖

 

屈曲位置為 PS 面靠近葉根整體葉片長度 5.2%處，由壓應力主導致使該位置出現(xiàn)局部屈曲進而導致整體屈曲。

 

我們選取葉片處展向長度 15.2%處及展向 30.4%處繪制其載荷位移曲線，如圖 4-5 所示，在OA 及 OC 段載荷系數較低，葉片各截面揮舞方向位移與加載力成線性關系，當載荷加載到達 A點及 C 點之后，位移對于加載力的響應更加敏感且為非線性關系，這表明在 A 點的載荷加載歷程中葉片出現(xiàn)了局部屈曲，當載荷繼續(xù)加載到 B、D 點后，在加載力不變的情況下?lián)]舞方向位移也會增加，Abaqus 停止繼續(xù)計算，表明此時已達到葉片的整體屈曲。

 

 

圖 4-5 葉片展向長度 15.2%及 30.4%位置加載力-位移關系圖

 

為更直觀的討論葉片局部屈曲與整體屈曲，比較局部屈曲位置單元結點與結構穩(wěn)定區(qū)域單元節(jié)點(隨機選取展向長度 54.3%位置主梁區(qū)域的節(jié)點)的應變及應力演化歷程，如圖 4-6 至圖 4-9。

 

 

圖 4-6 屈曲位置應變-弧長曲線

 

 

圖 4-7 葉片展向長度 54.3%處應變-弧長曲線

 

 

圖 4-8 屈曲位置應力-弧長曲線

 

 

圖 4-7 葉片展向長度 54.3%處應力-弧長曲線

 

通過上述對比可以發(fā)現(xiàn)，屈曲位置節(jié)點在弧長為1.2左右應力應變均有明顯的突變，而穩(wěn)定部位節(jié)點雖有斜率的增加，但變化微小。此外，從圖 4-6 及 4-8 可以看出，葉根屈曲位置受兩個方向的壓應力，以及較高的剪切力，足以將單元壓潰，致使發(fā)生屈曲。

 

對比特征值屈曲以及 Riks 非線性方法，計算預測到的屈曲位置有顯著差異，且 Riks 法得到的載荷系數也較特征值法小，這主要是因為改進的 Riks 方法因考慮結構的幾何缺陷以及結構的非線性關系 ，得到的結果更加準確可靠，同時也表明特征值只提供可能發(fā)生屈曲的位置，對于大變形的計算，該方法適用性不理想。。

 

結論

基于 Abaqus有限元軟件采用線性特征值以及非線性弧長法對時代新材 90m級葉片揮舞工況進行了屈曲分析，得到的結論如下：

 

1. 通過特征值屈曲方法得到的屈曲因子為 1.978，滿足了 GL2015 及 IEC61400-5 標準的要求，驗證了本司葉片結構設計的可靠性。

 

2.通過非線性 Riks 方法考慮了肯能存在的制造缺陷以及結構的非線性，得到的屈曲因子為

1.653，同樣滿足葉片設計標準。

 

3. 通過線性及非線性兩種方法計算結果的對比，表明特征值只提供可能發(fā)生屈曲的位置，對于大變形、大撓度的計算，非線性屈曲更加準確。

 

4.應用 Riks 方法可較準確的捕捉到葉片出現(xiàn)局部屈曲的位置，對葉片結構設計具有重要指導意義，可在葉片結構校核工作中進行推廣。

 

資料來源：達索官方