風(fēng)電機組葉片覆冰狀態(tài)與動力學(xué)特性的定量關(guān)系研究*

文 | 陳志剛，李錄平，楊波，劉勝先，封江，龔妙

我國的疆域遼闊，氣象條件復(fù)雜，在北部地區(qū)的陸地和海上，南部地區(qū)的高海拔山區(qū)等，冬季的風(fēng)雪霜凍給風(fēng)電機組的運行帶來了極大危害。運行在冬季濕冷地區(qū)的風(fēng)電機組，在冬季經(jīng)常會出現(xiàn)嚴(yán)重覆冰現(xiàn)象。工作中的風(fēng)電機組葉片覆冰后，會破壞風(fēng)電機組葉片的動力特性與氣動特性，冰的質(zhì)量、剛度等性質(zhì)會影響葉片的材料特性，造成葉片過載及載荷不均勻，使其特性與設(shè)計初衷不符，造成風(fēng)電機組不能正常發(fā)電，甚至導(dǎo)致嚴(yán)重的安全事故。及時、準(zhǔn)確地檢測出葉片的覆冰狀態(tài)，便于對風(fēng)電機組采取及時有效的運行和維護措施，對提高風(fēng)電場運行的安全可靠性和經(jīng)濟性，具有重要意義。

自20世紀(jì)90年代開始，對于構(gòu)件覆冰檢測的研究已經(jīng)引起人們廣泛關(guān)注，研究對象主要集中在航空領(lǐng)域（如飛機），主要的測量方式包括：機械法、電學(xué)法、熱學(xué)法、光學(xué)法、波導(dǎo)法等。但是，目前直接用于風(fēng)電機組葉片覆冰檢測的研究仍然比較少，能夠成熟用于工程實際的風(fēng)電機組葉片覆冰檢測方法與技術(shù)仍少有資料報道。

本文探索一種利用葉片模態(tài)參數(shù)變化特征來檢測葉片覆冰狀態(tài)的方法：以翼展1120mm（不包含葉片根部）長的小風(fēng)電機組葉片為研究對象，在Solidworks平臺上建立起三維實體模型。建立在振動力學(xué)、動力學(xué)分析以及柔性多體動力學(xué)理論基礎(chǔ)上，在ANSYS結(jié)構(gòu)分析軟件中對該葉片在覆冰前后狀態(tài)下模態(tài)進行分析，得到前6階固有頻率和振型，比較覆冰前后風(fēng)電機組葉片動力特性的變化情況，提取葉片覆冰后的振型曲率變化特征，獲得診斷葉片覆冰程度、覆冰位置的基本方法，為工程實際中診斷葉片覆冰故障提供理論依據(jù)。

研究對象簡介

本文以一小型風(fēng)電機組葉片為研究對象，在三維制圖工具軟件Solidworks中建立風(fēng)電機組葉片三維實體模型（見圖1）。如圖1所示，葉片為弦長140mm的等截面直葉片，內(nèi)部為實體結(jié)構(gòu)，長度為1120mm，由一種翼型界面拉伸而成。

圖1 模擬風(fēng)電機組葉片三維實體模型

圖2 葉片結(jié)冰沿翼展方向分布示意

葉片覆冰狀態(tài)描述模型

風(fēng)電機組葉片覆冰的空間分布描述采用如圖2所示的模型，假設(shè)葉片在某一區(qū)段附著一定厚度的冰，其覆冰狀態(tài)（Blade Icing Condition，簡稱BIC）信息用下式描述：

式中，x是葉片覆冰段沿翼展方向的起始坐標(biāo)（相對坐標(biāo)值）；Δx是葉片覆冰段沿翼展方向的長度（相對值）；λ是葉片覆冰段的覆冰厚度（相對坐標(biāo)值）。

一、覆冰區(qū)域起始點坐標(biāo)參數(shù)

定義葉片覆冰的起點坐標(biāo)參數(shù)為x（見圖2），其表示的是沿翼展方向覆冰區(qū)域的起始點相對位置，參數(shù)為一個無量綱量。以葉片與輪轂的結(jié)合處為翼展方向坐標(biāo)原點，定義描述覆冰起始位置的坐標(biāo)參數(shù)計算公式為：

式中，L是所測風(fēng)電機組葉片的翼展長度，m；X是覆冰段起始點距離葉根的坐標(biāo)，m；x的取值范圍為：0～1。

二、覆冰區(qū)域長度參數(shù)

定義覆冰長度參數(shù)為Δx，其表示在翼展方向覆冰區(qū)域沿翼展方向的相對長度。設(shè)覆冰段起點端到另一端的長度為 ΔX（見圖2），定義描述覆冰區(qū)域相對長度參數(shù)的計算公式為：

式中，ΔX是覆冰段沿翼展方向的長度，m。Δx的取值范圍為：0～1。

三、覆冰厚度參數(shù)

根據(jù)風(fēng)電場的風(fēng)能資源條件，風(fēng)電機組葉片的結(jié)構(gòu)大小多種多樣，因此用冰層厚度的絕對值不能準(zhǔn)確反映風(fēng)電機組葉片覆冰厚度。因此，定義覆冰厚度相對值λ，用其描述覆冰區(qū)域平均覆冰厚度：

式中，d是所測風(fēng)電機組葉片覆冰區(qū)域平均覆冰厚度，m。

覆冰葉片動力學(xué)特性計算數(shù)學(xué)模型

一、覆冰葉片等效材料特性參數(shù)計算模型

采用有限元方法計算覆冰葉片的動力學(xué)特性前，需要對比覆冰前與覆冰后葉片材料性能參數(shù)的變化情況。為此，本文在計算葉片中取中間一段180mm長的葉片段，建立葉片片段覆冰模型，覆冰相對厚度λ取值范圍為0～8.035714（1～9mm厚度）。根據(jù)所設(shè)置的葉片材料密度和冰層密度，用Solidworks軟件計算出葉片結(jié)構(gòu)的質(zhì)量變化情況。圖3、圖4所示為葉片段相對覆冰厚度λ分別取值為0.0、2.0時的覆冰模型。

圖3 葉片片段λ為0的模型

圖4 葉片片段λ為2的覆冰模型

圖5 葉片片段結(jié)構(gòu)網(wǎng)格化圖

圖6 葉片片段λ為2的結(jié)構(gòu)網(wǎng)格化圖

將覆冰前后的葉片段結(jié)構(gòu)導(dǎo)入ANSYS中，分別進行力學(xué)分析，在覆冰葉片段組合體中定義葉片材料的彈性模量E＝4.26e10，泊松比為0.22，密度為1950kg/m3；冰的彈性模量E＝2e9，泊松比為0.1，密度為900kg/m3。由于計算對象為小葉片，所以定義葉片為Solid45實體模型，網(wǎng)格化以后得到圖 5－圖6。

定義葉片段的一端為全約束，在葉片的另一端Y軸反方向施加一個50N力載荷，載荷施加圖和軟件計算后的結(jié)構(gòu)變形圖如圖7、圖8所示。

在葉片段的一端設(shè)置全約束，另一端施加一個力偶載荷產(chǎn)生5N?m力矩，通過ANSYS軟件計算得到覆冰葉片結(jié)構(gòu)的扭轉(zhuǎn)變形圖，如圖9所示。通過施加力和扭矩（力偶）分別得到葉片的彎曲變形和扭轉(zhuǎn)變形，在施加外界載荷中力和扭矩不變的條件下，分別對不同覆冰相對厚度的葉片進行結(jié)構(gòu)受力變形計算。

本文采用的抗彎剛度實驗測量公式為：

采用的扭轉(zhuǎn)剛度實驗測量公式為：

式中，F(xiàn)為結(jié)構(gòu)所受力載荷，N；M為結(jié)構(gòu)所受力偶載荷， N?m；y為結(jié)構(gòu)在力載荷下的位移，m；θ為結(jié)構(gòu)在力偶載荷下的扭轉(zhuǎn)角度，deg。

首先根據(jù)式（5）、（6）計算得到覆冰葉片片段結(jié)構(gòu)的抗彎剛度和扭轉(zhuǎn)剛度，結(jié)合質(zhì)量變化、抗彎剛度及扭轉(zhuǎn)剛度計算得到三項參數(shù)變化率隨覆冰厚度增加的變化關(guān)系圖，如圖10所示。

從圖10中可以看出，葉片的相對覆冰厚度λ從0到8之間變化時，葉片材料等效性能參數(shù)呈現(xiàn)如下變化規(guī)律：

（1）葉片覆冰后，材料等效質(zhì)量、等效彎曲剛度、等效扭轉(zhuǎn)剛度均有不同程度增加；

（2）葉片的等效質(zhì)量隨覆冰的相對厚度幾乎呈線性增加，當(dāng)相對覆冰厚度λ達(dá)到8時，葉片的等效質(zhì)量增加近35%左右；

（3）當(dāng)覆冰相對厚度λ在 4.5以內(nèi)時，葉片的等效彎曲剛度隨覆冰厚度的增加呈現(xiàn)快速增加趨勢，之后，隨著覆冰相對厚度的繼續(xù)增加，等效彎曲剛度增加的趨勢減緩，當(dāng)覆冰厚度達(dá)到8時，葉片的等效彎曲剛度增加近35%左右；

（4）當(dāng)覆冰相對厚度λ在4.5以內(nèi)時，葉片的等效扭轉(zhuǎn)剛度隨覆冰厚度的增加呈現(xiàn)緩慢增加趨勢，之后，隨著覆冰相對厚度的繼續(xù)增加，等效扭轉(zhuǎn)剛度增加的趨勢加快，當(dāng)覆冰厚度λ達(dá)到8.0時，葉片的等效扭轉(zhuǎn)剛度增加近12%左右。

圖7 葉片片段λ為0的彎曲變形圖

圖8 葉片片段λ為2的彎曲變形圖

圖9 葉片片段λ為2的扭轉(zhuǎn)變形圖

圖10 葉片片段在不同覆冰相對厚度下的材料性能變化圖

二、覆冰葉片模態(tài)參數(shù)變化量計算模型

（一）葉片模態(tài)頻率變化量計算模型

將葉片在揮舞方向的振動問題簡化成在根部固接的懸臂彎曲振動問題。懸臂梁的模態(tài)頻率計算問題可演化成由下式描述的特征值問題：

式中，K和M分別是葉片整體剛度矩陣和質(zhì)量矩陣；φ是正則化振型；ω是模態(tài)頻率。假定風(fēng)電機組葉片覆冰使葉片結(jié)構(gòu)剛度矩陣和質(zhì)量矩陣都產(chǎn)生一個細(xì)微攝動量，則對φ和ω也會產(chǎn)生一個細(xì)微改變量，改變后的振動方程為：

式中，ΔM，ΔK和Δφ分別為葉片結(jié)構(gòu)的質(zhì)量矩陣、剛度矩陣和振型的改變量。

對風(fēng)電機組葉片這樣的變截面懸臂梁結(jié)構(gòu)，在不同的覆冰厚度情況下，覆冰對質(zhì)量矩陣的影響程度和剛度矩陣影響程度是不同的。對式（8）適當(dāng)變形，得到：

在某一階振型下，有：

式（9）、（10）中，φT表示正則化振型矩陣φ的轉(zhuǎn)置矩陣。從上兩式可以分析得出，若剛度增加則葉片結(jié)構(gòu)自振頻率增加，剛度降低則結(jié)構(gòu)頻率會降低；而質(zhì)量增加，結(jié)構(gòu)頻率就會降低，反之，結(jié)構(gòu)頻率會增加。

（二）葉片模態(tài)振型曲率計算模型

由于風(fēng)電機組葉片為懸臂梁結(jié)構(gòu)，因為揮舞方向的抗彎能力最差，所以變形主要以揮舞方向彎曲變形為主，因此采用梁結(jié)構(gòu)單元進行推導(dǎo)，其振動微分方程為：

式中，x表示梁的軸向方向的坐標(biāo)；I（x）表示梁的抗彎截面模量；E、ρ分別表示材料的楊氏模量和密度；A（x）表示梁的截面積；C（x）表示材料的阻尼；y（x,t）表示梁的變形（位移）；f（x,t）表示作用在梁上的外力。

根據(jù)模態(tài)分析理論，式（11）的解可表示為各階模態(tài)振型的疊加形式：

式中，φi（x）為位移模態(tài)振型；qi（t）為模態(tài)坐標(biāo)；i為振型階數(shù)。

根據(jù)梁結(jié)構(gòu)單元振動微分方程式（11）和振動位移公

式（12），將式（11）化簡得到：

式中，ms、cs、ks分別表示梁的第s階模態(tài)質(zhì)量、阻尼和剛度；φs表示梁的第s階模態(tài)振型；qs表示梁的第s階模態(tài)坐標(biāo)。

根據(jù)材料力學(xué)給出的梁結(jié)構(gòu)靜力彎曲關(guān)系：

式中：M是葉片彎矩，ρ是曲率半徑。

設(shè)f（x,t）＝F（x）ejωt，由上兩式得到 ：

曲率是位移的二階導(dǎo)數(shù)，將上式離散化后得到：

矩陣形式為：

其中：

y′′的微分增量為 ：

覆冰葉片模態(tài)特性變化規(guī)律計算與分析

一、覆冰葉片模態(tài)頻率變化量計算與分析

在Solidworks中建立葉片覆冰模型。圖11－圖13分別表示在葉片不同覆冰起點位置上覆上一定長度區(qū)域的冰，并對不同覆冰狀態(tài)下的葉片進行模態(tài)分析。

采用ANSYS軟件，對不同覆冰情況下的葉片進行模態(tài)分析，得到頻率以及振型數(shù)據(jù)。葉片與冰的材料參數(shù)定義同上，將所有頻率數(shù)據(jù)整理出來，計算出相對于未覆冰狀態(tài)下葉片結(jié)構(gòu)的頻率變化率（圖14）。圖14中曲線說明有三個參數(shù)，第一個參數(shù)表示覆冰起點x，第二個參數(shù)表示覆冰長度Δx，第三個參數(shù)表示相對覆冰厚度λ。

圖14中包含所測葉片結(jié)構(gòu)前六階模態(tài)頻率。其中，從第一階頻率至第六階頻率分別為揮舞一階振型頻率、擺振一階振型頻率、揮舞二階振型頻率、揮舞三階振型頻率、扭轉(zhuǎn)一階振型頻率和揮舞四階振型頻率。

圖11 覆冰起點x＝0.25，覆冰長度Δx＝0.5，覆冰厚度λ＝4時葉片三維圖

圖12 覆冰起點x＝0.25，覆冰長度Δx＝0.25，覆冰厚度λ＝4時葉片三維圖

圖13 覆冰起點x＝0.75，覆冰長度Δx＝0.25，覆冰厚度λ＝2時葉片三維圖

從圖14可以發(fā)現(xiàn)：

（1）葉片表面覆冰后，會導(dǎo)致模態(tài)頻率的變化；覆冰越嚴(yán)重，模態(tài)頻率的相對變化值越大。

（2）不同覆冰狀態(tài)下各階模態(tài)頻率的變化率也是不同的。整體而言，葉片覆冰對一階模態(tài)頻率影響最大。

（3）當(dāng)葉片所有表面覆冰時，葉片模態(tài)頻率變化最大且與覆冰厚度成正比關(guān)系。

（4）葉片根部覆冰后導(dǎo)致葉片結(jié)構(gòu)剛度增加明顯，故模態(tài)頻率升高；葉片葉尖區(qū)域覆冰后導(dǎo)致葉片結(jié)構(gòu)質(zhì)量增加明顯，引起模態(tài)頻率下降。

圖14 不同覆冰狀態(tài)下葉片結(jié)構(gòu)前六階模態(tài)頻率變化

圖15 葉片各點在不同覆冰條件下一階曲率變化率

二、葉片覆冰狀態(tài)與葉片模態(tài)振型的關(guān)系

由于葉片覆冰后，結(jié)構(gòu)剛度顯著增大，導(dǎo)致的曲率變化［Δy′′］，主要由［Δφ′′r］［ΔYr］［Δφr］三者共同作用。而［Δy′′］ 與［Δφ′′r］的變化在葉片位置坐標(biāo)上存在統(tǒng)一對應(yīng)性的變化。

用ANSYS軟件仿真計算出的葉片結(jié)構(gòu)固有特性參數(shù)包括模態(tài)頻率和模態(tài)振型，所獲得的模態(tài)振型主要以位移模態(tài)數(shù)據(jù)為主。根據(jù)差分法計算出葉片結(jié)構(gòu)的各階曲率模態(tài)參數(shù)。這里主要分析葉片結(jié)構(gòu)的揮舞方向的模態(tài)曲率，根據(jù)下式得出葉片結(jié)構(gòu)的曲率變化率曲線圖，如圖15所示。

式中：ξi為i點處覆冰后的曲率變化率；φwi為i點處覆冰后曲率值；φfi為i點處無冰狀態(tài)下曲率值。

從圖15中可以看到：

（1）在葉片全覆冰情況下，覆冰2mm時曲率變化較小，曲率變化率曲線無規(guī)則；當(dāng)覆冰達(dá)到4mm時曲率變化明顯增大，同時曲率變化率曲線呈現(xiàn)明顯的拋物線形狀。由此可見，在覆冰區(qū)域一定時，覆冰厚度越大，葉片模態(tài)振型曲率的變化率越大。

（2）覆冰起點為0.5，覆冰長度分別為0.25和0.50時，葉片對應(yīng)位置的曲率變化率曲線也有凸起形狀，同時可以看出覆冰長度Δx＝0.50的曲率變化率大于覆冰長度Δx＝0.25的曲率變化率。由此可見，在覆冰起點位置相同時，覆冰區(qū)域越大，葉片模態(tài)振型曲率的變化率越大。

從上述分析結(jié)論可以看出，通過對覆冰前后葉片結(jié)構(gòu)的曲率變化進行分析研究，能夠進一步對覆冰狀態(tài)進行判斷。

攝影：李博

結(jié)論

本文提出了風(fēng)電機組葉片覆冰狀態(tài)的空間描述模型，給出了葉片覆冰起點參數(shù)、覆冰長度參數(shù)和覆冰厚度參數(shù)的定義；給出了覆冰葉片動力學(xué)參數(shù)變化量計算模型。以一小型風(fēng)電機組葉片為例，在Solidworks三維建模軟件中建立三維實體模型，應(yīng)用ANSYS結(jié)構(gòu)分析軟件對風(fēng)電機組葉片進行有限元模型分析，通過有限元計算，獲得葉片在不同覆冰狀態(tài)下的模態(tài)頻率變化率、模態(tài)振型曲率變化率，通過分析發(fā)現(xiàn)：

（1）葉片覆冰區(qū)域靠近葉尖部分時，葉片模態(tài)頻率呈現(xiàn)下降趨勢；葉片所有表面覆冰時，模態(tài)頻率變化率與覆冰厚度成近似正比關(guān)系。葉片根部區(qū)域覆冰后，葉片結(jié)構(gòu)以剛度增加為主，故模態(tài)頻率升高。葉片葉尖區(qū)域覆冰后，葉片結(jié)構(gòu)以質(zhì)量分布變化為主，導(dǎo)致模態(tài)頻率下降。

（2）葉片模態(tài)振型曲率變化率值與葉片覆冰位置之間關(guān)系密切，根據(jù)葉片模態(tài)振型曲率變化率的不同，能夠在一定程度上對不同覆冰情況下葉片覆冰位置和覆冰厚度進行判斷。

（3）葉片模態(tài)頻率變化率、模態(tài)振型曲率變化率與葉片覆冰狀態(tài)具有明顯的定量關(guān)系。工程實際中，通過這種定量關(guān)系，可以診斷葉片的覆冰故障。