武浩 來永斌 王龍 李亮

摘 要：針對風(fēng)力機葉片在正常工況下運行時受到周期性的氣動力導(dǎo)致葉片發(fā)生振動，降低葉片使用壽命的情況，研究了風(fēng)力機葉片在不同風(fēng)速下的振動特性。選取不同風(fēng)速條件下的5種工況 （風(fēng)速范圍為15～40 m/s），選用CFD方法對NREL PHASE VI葉片進(jìn)行模擬計算，獲取不同風(fēng)速下的振型和振動位移曲線。結(jié)果表明：葉片的主要振型是揮舞和擺振，高階葉片振型存在著彎曲和扭轉(zhuǎn)組合的復(fù)雜變形；來流速度從15 m/s增大到40 m/s時，葉片吸力面的壓力分布不均勻性不斷提高，來流速度為40 m/s時最大壓力差值約達(dá)到3 000 Pa；來流速度為15 m/s時振幅最小為0.525 4 mm，來流速度為40 m/s時振幅最大，為3.628 2 mm，約是最小振幅的6.9倍；5種工況的振動曲線均呈現(xiàn)衰減趨勢，葉片趨于穩(wěn)定振動；當(dāng)來流風(fēng)速越大時，由來流風(fēng)所產(chǎn)生的氣動力對葉片的作用力越大，葉片的振幅呈現(xiàn)增大的趨勢。研究結(jié)果可為風(fēng)力機設(shè)計提供參考。

關(guān)鍵詞：非線性振動力學(xué)；非線性振動系統(tǒng)；風(fēng)力機葉片；計算流體力學(xué)；CFD

中圖分類號：TK83 文獻(xiàn)標(biāo)志碼：A

文章編號：1008-1542（2018）05-0401-08

風(fēng)力機葉片是風(fēng)電機組獲取風(fēng)能的關(guān)鍵部件，也是承受載荷復(fù)雜的部件之一，其氣動性能決定了風(fēng)力機的風(fēng)能轉(zhuǎn)化率和安全運行穩(wěn)定性。在不同來流風(fēng)速下，風(fēng)力機葉片會產(chǎn)生振動，葉片會不斷從氣流吸收能量，當(dāng)此能量大于結(jié)構(gòu)阻尼負(fù)功時，葉片振動會不斷加劇，產(chǎn)生顫振[1-2]，導(dǎo)致風(fēng)力機損壞。因此，研究風(fēng)力機葉片在不同來流速度下的振動特性，對于風(fēng)力機安全運行具有重要意義。

國內(nèi)外學(xué)者對風(fēng)力機葉片氣彈穩(wěn)定問題進(jìn)行了相關(guān)研究。KOTTAPLLI等[3]用非耦合的旋轉(zhuǎn)模態(tài)對葉片氣動穩(wěn)定性進(jìn)行了研究，得出了風(fēng)力機葉片的靜態(tài)響應(yīng)及其穩(wěn)定性邊界。KIM等[4]采用非線性有限元方法，在考慮旋轉(zhuǎn)效應(yīng)的基礎(chǔ)上，對5 MW風(fēng)力機葉片進(jìn)行了氣動性能預(yù)測。BAZILEVS等[5]和HSU等[6]對NREL 5 MW海上風(fēng)力機進(jìn)行雙向流固耦合數(shù)值計算研究，有利于海上風(fēng)力機葉片的結(jié)構(gòu)設(shè)計優(yōu)化和新型材料的運用。張瑞琴等[7]使用ANSYS分析NACA0012翼型在流固耦合作用下的顫振特性，結(jié)果表明來流速度和攻角是影響葉片氣動彈性穩(wěn)定性的重要因素。WANG等[8]基于有限體積方法與LU-SGS算法相結(jié)合以分析風(fēng)力機葉片振動特性。李亮等[9]分析風(fēng)力機葉片揮舞-擺振氣彈穩(wěn)定性，結(jié)果表明通過設(shè)置安裝角，利用揮舞-擺振耦合效應(yīng)可以改善葉片氣彈穩(wěn)定性。任勇生等[10]采用特征值方法和時域積分法估計葉片的顫振邊界，對NACA0012風(fēng)力機葉片進(jìn)行數(shù)值分析，結(jié)果表明風(fēng)力機在高風(fēng)速運行時，葉片作發(fā)散運動且系統(tǒng)是不穩(wěn)定的。李本立等[11]建立了風(fēng)力機轉(zhuǎn)子葉片的非線性運動方程，采用模態(tài)法求解揮舞、擺振、扭轉(zhuǎn)微分方程，利用數(shù)值結(jié)果對風(fēng)力機的氣動彈性穩(wěn)定性進(jìn)行了分析。杜標(biāo)等[12]利用S-A湍流模型對S809翼型進(jìn)行動態(tài)失速分析，結(jié)果表明升力系數(shù)與風(fēng)速呈線性增長關(guān)系。葉學(xué)民等[13]利用雙凹槽葉頂結(jié)構(gòu)，模擬風(fēng)機在不同開槽深度下的性能，結(jié)果表明葉頂開槽深度為3 mm時效果最好?？芎＼姷萚14]通過振動理論和有限元分析對航空發(fā)動機葉片進(jìn)行靜頻、動頻和振動響應(yīng)分析諧響應(yīng)分析，結(jié)果表明在葉片前緣最容易發(fā)生破壞。李克安等[15]利用Ritz-Galerkin方法對航空發(fā)動機轉(zhuǎn)子葉片靜頻和動頻進(jìn)行計算分析，結(jié)果表明該方法較適合流固耦合振動分析。王立剛等[16]依據(jù)轉(zhuǎn)子理論，結(jié)合Runge-Kutta研究葉片振動對轉(zhuǎn)子-軸承系統(tǒng)動力學(xué)行為的影響，結(jié)果表明葉片結(jié)構(gòu)阻尼對該系統(tǒng)動力學(xué)行為影響較大。通過上述分析可知，對于高風(fēng)速下風(fēng)力機葉片振動的研究較少。

筆者采用ANSYS Workbench對5種不同風(fēng)速條件下（風(fēng)速范圍為15～40 m/s）的NREL PHASE VI葉片進(jìn)行流固耦合仿真，分析了風(fēng)力機葉片的振型和振動特性。

1 控制方程及湍流模型

1.1 結(jié)構(gòu)動力學(xué)方程

2 結(jié)構(gòu)模型、參數(shù)設(shè)置及網(wǎng)格劃分

2.1 計算模型

研究對象為NREL PHASE VI葉片。材料屬性設(shè)置：彈性模量[WTBX]E=73 GPa，密度ρ=2 540 kg/m3，泊松比λ=0.22。葉片模型如圖1所示。流體為標(biāo)準(zhǔn)狀態(tài)下的空氣。對葉片實體進(jìn)行網(wǎng)格劃分，共劃分2 800個六面體網(wǎng)格。

2.2 模擬參數(shù)設(shè)置及網(wǎng)格劃分

由于其工作環(huán)境為近海區(qū)域，所以模擬單個葉片不轉(zhuǎn)動的情況。工作環(huán)境的來流速度分別為15，20，25，30和40 m/s等5種工況。圖2為流場計算區(qū)域網(wǎng)格。整個計算域均采用四面體非結(jié)構(gòu)化網(wǎng)格，共計約642 000個網(wǎng)格單元。計算域在進(jìn)口邊界設(shè)置不同來流速度，其來流風(fēng)湍流分?jǐn)?shù)強度（fractional intensity）為5%，湍流長度尺寸（eddy length size）為0.25 m，分別對不同工況進(jìn)行計算分析。進(jìn)口采用速度進(jìn)口，出口采用opening邊界條件。

3 計算結(jié)果與分析

3.1 經(jīng)典算例驗證

為驗證所提出的流固耦合方法的可行性，用圓柱繞流振動的算例進(jìn)行驗證。

3.1.1 模型參數(shù)和邊界條件

模型結(jié)構(gòu)圖如圖3所示，流體域采用四面體非結(jié)構(gòu)網(wǎng)格，網(wǎng)格數(shù)共14萬，圓柱結(jié)構(gòu)有限元網(wǎng)格數(shù)為3 000。其中入口采用速度進(jìn)口；出口采用opening邊界條件。

3.1.2 結(jié)果驗證

圖4給出了壓力系數(shù)對比圖，其中虛線為本文算例所得壓力系數(shù)分布曲線，實線為相同條件下文獻(xiàn)[20]中的參考曲線。從圖4可以發(fā)現(xiàn)，兩條曲線的變化趨勢及數(shù)值大小基本吻合，且分離點分布均在80°左右，表明本文所采用的流固耦合計算方法是相對可行有效的。

3.2 葉片模態(tài)分析

固體計算選用ANSYS Workbench 瞬態(tài)動力學(xué)結(jié)構(gòu)模塊進(jìn)行分析，結(jié)構(gòu)瞬態(tài)分析和流場的非定常計算采用相同的時間步長，設(shè)置耦合步長為2×10-4 s，并施加約束條件。表1和圖5分別是葉片的前四階固有頻率及振型。

根據(jù)圖5可知，葉片振型主要有3種：揮舞、擺振和扭轉(zhuǎn)振動。其中揮舞是切向的彎曲振動；擺振是軸向的彎曲振動；扭振是徑向的扭轉(zhuǎn)振動。葉片的一階振型是[WTBX]Y向的揮舞振動，葉片的二階振型是擺振，葉片三階振型是Y向彎曲振動但相對于一階振型的彎曲幅度更大，葉片四[WT]階振型是彎曲和扭轉(zhuǎn)組合振動，振型較為復(fù)雜。根據(jù)結(jié)構(gòu)振動理論[20]，振動發(fā)生的主要能量集中在一階和二階，因此葉片的主要振型是揮舞和擺振，高階葉片振型存在著彎曲和扭轉(zhuǎn)組合的復(fù)雜變形。

3.3 葉片吸力面壓力分析

圖6顯示不同來流速度下葉片吸力面的壓力云圖。由圖6可知：當(dāng)來流速度為15 m/s時，除了葉根和葉中小部分位置，其余位置沿著葉展方向壓力呈現(xiàn)平行分布趨勢。伴隨著來流速度增大到20 m/s，葉片前緣的吸力峰值隨之變大，而且靠近葉根和葉中部分區(qū)域壓力變小，等壓值構(gòu)成的平行分布區(qū)域開始縮小，同時靠近葉尖位置處帶狀區(qū)域范圍也在縮小。從圖6 c）、d）可以看出，當(dāng)來流速度達(dá)到25 m/s和30 m/s時，壓力分布的不均勻性進(jìn)一步提高，靠近葉根和葉中區(qū)域的壓力更小，葉尖的葉片前緣處等壓區(qū)域壓力值增大。當(dāng)風(fēng)速達(dá)到40 m/s時，葉片吸力面前緣壓強達(dá)到最大值，壓力分布不均勻性達(dá)到最大程度，其中最大壓力差值約達(dá)到3 000 Pa。

3.4 不同來流速度對葉片振動的影響

本文選擇風(fēng)力機葉片工作環(huán)境的來流速度是15～40 m/s。圖7是來流速度為15 m/s時[WTBX]的Y向位移云圖，[WT]可以看出振動位移量主要產(chǎn)生在葉尖部分。根據(jù)圖8葉片壓力云圖可知，葉片的壓力主要集中在正對來流風(fēng)速的葉片前緣。在葉片葉尖前緣處壓力大且振動位移較大，因此振動曲線監(jiān)控采集點選取葉片前緣頂部。

根據(jù)表2和圖9可以看出，當(dāng)來流風(fēng)速越大時，由來流風(fēng)所產(chǎn)生的氣動力對葉片的作用力越大，葉片的振幅呈現(xiàn)增大的趨勢。其中來流速度為15 m/s時振幅最小為0.525 4 mm，來流速度為40 m/s時振動達(dá)到最大為3.628 2 mm，約是最小振幅的6.9倍；從30 m/s到40 m/s振幅增幅最大。

風(fēng)力機葉片在各種工況工作過程中，除了葉片本身所受的重力載荷和葉片工作過程所受的離心力之外，受到最主要的外載荷來自來流風(fēng)產(chǎn)生的氣動力。根據(jù)圖10 a）—e）可以看出5種工況的振動曲線呈現(xiàn)衰減的趨勢，葉片趨于穩(wěn)定振動。分析認(rèn)為，流體產(chǎn)生的氣動力作用于葉片引起葉片的變形和振動，同時由于葉片翼型會存在升力[21]，因此葉片會偏離初始位置先向上振動，當(dāng)葉片達(dá)到最大振幅時，受重力載荷和材料彈性力的作用，葉片開始向相反方向振動，由于慣性力和重力的作用，葉片經(jīng)過平衡位置繼續(xù)振動達(dá)到反向最大振幅。由于葉片振動消耗外界能量（結(jié)構(gòu)阻尼所做負(fù)功大于氣動力所做正功，整體呈現(xiàn)衰減狀態(tài)），振動最大位移與振動初始位置之間距離較小，并形成周期性衰減振動過程。

4 結(jié) 論

建立了NREL PHASE VI葉片在不同來流的流固耦合數(shù)值模型，得到了前4階固有頻率、振型和不同工況下的振動特征曲線。葉片的主要振型是揮舞和擺振，高階葉片振型存在著彎曲和扭轉(zhuǎn)組合的復(fù)雜變形。來流速度為15 m/s時振幅最小為0.525 4 mm，來流速度為40 m/s時振動達(dá)到最大為3.628 2 mm，約是最小振幅的6.9倍；5種工況的振動曲線呈現(xiàn)衰減的狀態(tài)，振動穩(wěn)定；當(dāng)來流風(fēng)速越大時，由來流風(fēng)所產(chǎn)生的氣動力對葉片的作用力越大，葉片的振幅呈現(xiàn)增大的趨勢。研究結(jié)果可為風(fēng)力機設(shè)計提供參考。本文所做工作集中在靜止?fàn)顟B(tài)的風(fēng)力機葉片，下一步將對不同轉(zhuǎn)速下的旋轉(zhuǎn)葉片開展相關(guān)振動研究，同時，擬擴大風(fēng)速范圍，以使研究結(jié)果具有更好的通用性。

參考文獻(xiàn)/References：

[1] 趙永輝. 氣動彈性力學(xué)與控制[M]. 北京：科學(xué)出版社， 2007.

[2] 周盛. 葉輪機氣動彈性力學(xué)引論[M].北京：國防工業(yè)出版社， 1989.

[3] KOTTAPALLI S，F(xiàn)RIEDMANN P， ROSEN A. Aeroelastic stability and response of horizontal axis wind turbine blades[J]. Aiaa Journal， 1979， 17（12）： 49-66.

[4] KIM D H， KIM Y H. Performance prediction of a 5 MW wind turbine blade considering aeroelastic effect[J]. Proceedings of World Academy of Science Engineering & Technology， 2013（81）：771-775.

[5] BAZILEVS Y， HSU M C， KIENDL J， et al. 3D simulation of wind turbine rotors at full scale（Part Ⅱ）： Fluid-structure interaction modeling with composite blades[J]. International Journal for Numerical Methods in Fluids， 2011， 65（1/2/3）：207-235.

[6] JONKMAN J， BUTTERFIELD S， MUSIAL W， et al. Definition of A 5-MW Reference Wind Turbine for Offshore System Development[R]. [S.l.]：Office of Scientific & Technical Information Technical Reports， 2009.

[7] 張瑞琴， 翁建生. 基于流固耦合的葉片顫振分析[J]. 計算機仿真， 2011， 28（3）： 48-51.

ZHANG Ruiqin， WENG Jiansheng. Analysis of blade flutter based on fluid solid coupling [J]. Computer Simulation， 2011， 28 （3）： 48-51.

[8] WANG Long， SUN Lunye， WU Guang， et al. Research on algorithm of blade vibration for general wind turbine[C]// International Symposium on Precision Mechanical Measurements. [S.l.]：[s.n.]， 2016：990325.

[9] 李亮， 李映輝， 楊鄂川. 風(fēng)力機葉片揮舞-擺振氣彈穩(wěn)定性分析[J]. 噪聲與振動控制， 2015， 35（5）： 30-34.

LI Liang， LI Yinghui， YANG Echuan. Aerodynamic stability analysis of wind turbine blades waving - shimmy [J]. Noise and Vibration Control， 2015， 35 （5）： 30-34.

[10]任勇生， 劉廷瑞， 楊樹蓮. 風(fēng)力機復(fù)合材料葉片的動力失速氣彈穩(wěn)定性研究[J]. 機械工程學(xué)報， 2011， 47（12）：113-125.

REN Yongsheng， LIU Tingrui， YANG Shulian. Study on aerodynamic stall aeroelastic stability of composite blades of wind turbine [J]. Journal of Mechanical Engineering， 2011， 47 （12）： 113-125.

[11]李本立， 安玉華. 風(fēng)力機氣動彈性穩(wěn)定性的研究[J]. 太陽能學(xué)報， 1996（4）：314-320.

LI Benli， AN Yuhua. Study on aeroelastic stability of wind turbines [J]. Acta Solar Energy， 1996（4）： 314-320.

[12]杜標(biāo)， 李雪斌， 王龍，等. 不同風(fēng)速下風(fēng)力機動態(tài)特性研究[J]. 制造業(yè)自動化， 2017，39（12）：86-94.

DU Biao， LI Xuebin， WANG Long， et al. Research on dynamic characteristics of wind turbines under different wind speeds [J]. Automation of Manufacturing Industry， 2017， 39（12）：86-94.

[13]葉學(xué)民， 李鵬敏， 李春曦. 雙凹槽葉頂結(jié)構(gòu)下的軸流風(fēng)機性能及葉片振動特性研究[J]. 機械工程學(xué)報， 2015， 51（4）：167-174.

YE Xuemin， LI Pengmin， LI Chunxi. Investigation of double grooved blade tip structure on aerodynamic performance and vibration characteristic of an axial flow fan[J]. Journal of Mechanical Engineering， 2015， 51（4）：167-174.

[14]寇海軍， 張俊紅， 林杰威. 航空發(fā)動機風(fēng)扇葉片振動特性分析[J]. 西安交通大學(xué)學(xué)報， 2014， 48（11）：109-114.

KOU Haijun， ZHANG Junlin， LIN Jiewei. Aero-engine fan blade vibration characteristic analysis[J]. Journal of

Xian Jiaotong Univer-sity， 2014， 48（11）：109-114.

[15]李克安， 林左鳴， 楊勝群，等. 航空發(fā)動機轉(zhuǎn)子葉片振動方程及其頻率計算[J]. 航空學(xué)報， 2013， 34（12）：2733-2739.

LI Kean， LIN Zuoming， YANG Shengqun， et al. Vibration equation and frequency calculation of aero engine rotor blades [J]. Acta Aeronautics Sinica， 2013， 34 （12）： 2733-2739.

[16]王立剛， 曹登慶， 胡超，等. 葉片振動對轉(zhuǎn)子-軸承系統(tǒng)動力學(xué)行為的影響[J]. 哈爾濱工程大學(xué)學(xué)報， 2007， 28（3）：320-325.

WANG Ligang， CAO Dengqing， HU Chao， et al. Effect of the blade vibration on the dynamical behaviors of a rotor-bearing system[J]. Journal of Harbin Engineering University， 2007， 28（3）： 320-325.

[17]胡海巖. 機械振動基礎(chǔ)[M]. 北京：北京航空航天大學(xué)出版社， 2005.

[18]周光垌. 流體力學(xué)（上冊）[M]. 北京：高等教育出版社， 1992.

[19]MENTER F R. Two-equation eddy-viscosity turbulence models for engineering applications[J]. Aiaa Journal， 1994， 32（8）：1598-1605.

[20]白樺， 李加武， 夏勇. 低雷諾數(shù)圓柱繞流數(shù)值模擬及控制措施[J]. 建筑科學(xué)與工程學(xué)報， 2010， 27（4）：39-43.

BAI Hua， LI Jiawu ， XIA Yong. Numerical simulation and control measures of flow around circular cylinders at low reynolds number[J]. Journal of Architecture & Civil Engineering， 2010， 27（4）： 39-43.

[21]翁智遠(yuǎn). 結(jié)構(gòu)振動理論[M]. 上海：同濟大學(xué)出版社， 1988.

[22]李建波， 高正. 涵道風(fēng)扇空氣動力學(xué)特性分析[J]. 南京航空航天大學(xué)學(xué)報， 2005， 37（6）：680-684.

LI Jianbo， GAO Zheng. Aerodynamic characteristics analysis of ducted fan[J]. Journal of Nanjing University of Aeronautics & Astronautics， 2005， 55（4）： 2373-2381.