小型垂直軸風(fēng)力機(jī)葉片流固耦合分析

韓兆輝，郭志平，張艷鋒

（1.內(nèi)蒙古工業(yè)大學(xué)機(jī)械工程學(xué)院，內(nèi)蒙古 呼和浩特 010051；2.內(nèi)蒙古工業(yè)大學(xué)理學(xué)院，內(nèi)蒙古 呼和浩特 010051）

1 引言

可再生能源在現(xiàn)今能源短缺的時(shí)代已經(jīng)充分顯示了本身所具有的潛力及重要性。風(fēng)能作為可再生能源的重要組成部分，使用風(fēng)力發(fā)電已經(jīng)成為世界上最具前景的能源開發(fā)方式。葉片作為風(fēng)力發(fā)電機(jī)的重要組成部件，直接影響風(fēng)力發(fā)電機(jī)的工作效率。在風(fēng)力機(jī)運(yùn)行過程中會(huì)引起葉片發(fā)生變形和振動(dòng)，隨著葉片變形和振動(dòng)的產(chǎn)生又會(huì)對(duì)流場(chǎng)中的氣流產(chǎn)生影響，為典型流固耦合問題，所以對(duì)葉片與周圍流場(chǎng)進(jìn)行耦合作用分析十分重要。

國外研究人員對(duì)風(fēng)力發(fā)電機(jī)組的研究早于國內(nèi)。文獻(xiàn)［1］對(duì)翼型的流固耦合理論做了系統(tǒng)總結(jié)；文獻(xiàn)［2］對(duì)偏航條件下的大型風(fēng)力機(jī)進(jìn)行了仿真計(jì)算，對(duì)該風(fēng)機(jī)的氣動(dòng)性能進(jìn)行了一系列的計(jì)算分析；文獻(xiàn)［3］基于流固耦合方法對(duì)風(fēng)力機(jī)葉片進(jìn)行了仿真分析，對(duì)單個(gè)葉片的應(yīng)變分布進(jìn)行了計(jì)算分析；文獻(xiàn)［4］對(duì)大型風(fēng)力機(jī)葉片進(jìn)行了計(jì)算分析，考慮氣動(dòng)彈性下葉片的影響，并對(duì)大型風(fēng)力機(jī)葉片性能進(jìn)行了理論分析。

80年代中期開始，國內(nèi)研究人員開始對(duì)風(fēng)力發(fā)電機(jī)流固耦合等問題進(jìn)行研究［5］。文獻(xiàn)［6］針對(duì)離網(wǎng)小型垂直軸風(fēng)力發(fā)電機(jī)，對(duì)風(fēng)力發(fā)電機(jī)葉片進(jìn)行了流固耦合仿真分析，研究了翼型的氣動(dòng)特性；文獻(xiàn)［7］對(duì)風(fēng)機(jī)建筑一體化中垂直軸風(fēng)力機(jī)葉片和主軸的受力情況，使用ANSYS Workbench運(yùn)用單項(xiàng)流固耦合的方法，在不同來流條件下對(duì)葉片和主軸的靜應(yīng)力進(jìn)行了分析和比較；文獻(xiàn)［8］針對(duì)在ANSYS Workbench軟件中不能在Fluent模塊下直接觀察到葉片在流場(chǎng)中所受的壓力，通過軟件ANSYS-CFX 組合，找到分析流固耦合的方法，使風(fēng)輪數(shù)值計(jì)算分析更合理、更準(zhǔn)確；文獻(xiàn)［9］引入類似結(jié)構(gòu)方程組處理流體區(qū)域連續(xù)變形，對(duì)流體網(wǎng)格節(jié)點(diǎn)坐標(biāo)進(jìn)行更新，采用強(qiáng)耦合方法對(duì)垂直軸風(fēng)力機(jī)葉片進(jìn)行了流固耦合分析，研究了葉片的氣動(dòng)性能。

流固耦合法是觀察分析固體場(chǎng)在流體場(chǎng)作用下的一系列行為及固體場(chǎng)反作用于流體場(chǎng)時(shí)流體場(chǎng)的變化情況。應(yīng)用ANSYS Workbench仿真軟件對(duì)垂直軸風(fēng)力機(jī)葉片進(jìn)行雙向流固耦合分析，是對(duì)仿真葉片進(jìn)行模擬分析的重要手段。所得結(jié)果可以為葉片優(yōu)化設(shè)計(jì)、研發(fā)以及安全運(yùn)轉(zhuǎn)提供一定的參考。

2 流固耦合理論基礎(chǔ)及方法

2.1 流體運(yùn)動(dòng)守恒方程

流體運(yùn)動(dòng)必須遵循基本的方程，包括連續(xù)性方程、動(dòng)量方程-牛頓運(yùn)動(dòng)定律、能量方程-熱力學(xué)第一定律。這三個(gè)守恒方程是CFD的理論基石。針對(duì)這里所涉及的問題，連續(xù)性方程和動(dòng)量方程可通過如下方程進(jìn)行描述。

連續(xù)性方程為：

式中：t—時(shí)間；ρ—流體密度；vx、vy、vz—沿X軸、Y軸、Z軸的速度。

任何流體問題都必須滿足質(zhì)量守恒定律，即單位時(shí)間內(nèi)流體微元體中質(zhì)量的增加，等于同一時(shí)間間隔內(nèi)流入該微元的凈質(zhì)量［10］。動(dòng)量方程-牛頓運(yùn)動(dòng)定律：

式中：F—合力；m—質(zhì)量；v—速度；ρ—密度；

2.2 流固耦合方程

在風(fēng)力機(jī)運(yùn)轉(zhuǎn)過程中葉片與流場(chǎng)之間相互作用。葉片受到來自流場(chǎng)壓力的同時(shí)會(huì)發(fā)生變形，變形后的葉片反作用于流場(chǎng)，進(jìn)一步擾動(dòng)葉片所處的流場(chǎng)。因此，在流固耦合壁面處應(yīng)滿足以下方程。

式中：τ—應(yīng)力；d—位移；q—熱流量；T—溫度。

3 三維模型的建立

這里利用SolidWorks建模軟件對(duì)風(fēng)力機(jī)葉片進(jìn)行三維實(shí)體模型的建立。風(fēng)力機(jī)葉片三維實(shí)體模型與室外實(shí)驗(yàn)場(chǎng)地使用的風(fēng)力機(jī)葉片為同一規(guī)格型號(hào)。這里所使用仿真葉片翼型是NACA0021，該仿真葉片內(nèi)部由三根鋁管加強(qiáng)筋貫穿五片支撐片，葉片兩端由擋板固定，外表進(jìn)行蒙皮，使用鉚釘固定組裝而成。建立的仿真葉片模型總長(zhǎng)h=1200mm，弦長(zhǎng)c=265mm，重量G=4.28kg，最后進(jìn)行葉片表面處理，這樣使仿真葉片模型趨近于室外實(shí)驗(yàn)所用模型，且對(duì)葉片強(qiáng)度以及剛度都有所增強(qiáng)。具體模型，如圖1所示。

圖1 葉片整體模型Fig.1 Blade Integral Model

4 實(shí)驗(yàn)設(shè)備、測(cè)試方案及結(jié)果

4.1 實(shí)驗(yàn)設(shè)備

該實(shí)驗(yàn)使用的風(fēng)洞隸屬于內(nèi)蒙古自治區(qū)風(fēng)能太陽能利用機(jī)理及優(yōu)化重點(diǎn)實(shí)驗(yàn)室［11］。實(shí)驗(yàn)對(duì)象為小型垂直軸風(fēng)力機(jī)葉片，實(shí)驗(yàn)采用位移測(cè)試分析系統(tǒng)對(duì)葉片進(jìn)行耦合實(shí)驗(yàn)，系統(tǒng)主要包括位移傳感器、位移信號(hào)采集分析軟件、加速度傳感器和數(shù)據(jù)線等配套設(shè)備組裝而成。

4.2 測(cè)試方案

利用位移測(cè)試分析系統(tǒng)在垂直軸風(fēng)力機(jī)研究室進(jìn)對(duì)實(shí)驗(yàn)葉片進(jìn)行位移測(cè)試，對(duì)實(shí)驗(yàn)葉片進(jìn)行監(jiān)測(cè)點(diǎn)布置，在葉片端部尾緣、0.25h尾緣與0.5h尾緣處分別選取一個(gè)監(jiān)測(cè)點(diǎn)，定義為監(jiān)測(cè)點(diǎn)一、二、三。在葉片監(jiān)測(cè)點(diǎn)處布置與實(shí)驗(yàn)葉片表面方向垂直的位移傳感器，對(duì)位移傳感器編號(hào)與其對(duì)應(yīng)的監(jiān)測(cè)點(diǎn)編號(hào)一致，選用實(shí)驗(yàn)風(fēng)速為8m/s對(duì)葉片進(jìn)行實(shí)驗(yàn)分析。實(shí)驗(yàn)所用支架為剛性結(jié)構(gòu)，支架與葉片、地面采用螺栓連接。葉片置于風(fēng)洞出風(fēng)口處中間位置，弦長(zhǎng)與風(fēng)速垂直。三維實(shí)驗(yàn)測(cè)試簡(jiǎn)圖，如圖2所示。

圖2 三維實(shí)驗(yàn)測(cè)試簡(jiǎn)圖Fig.2 Three Dimensional Experimental Test Diagram

4.3 實(shí)驗(yàn)結(jié)果

通過在不同工況下對(duì)葉片進(jìn)行位移實(shí)驗(yàn)，得到位移傳感器X、Y、Z3個(gè)方向上的位移信號(hào)，對(duì)信號(hào)放大后進(jìn)行分析處理，并對(duì)信號(hào)進(jìn)行一系列函數(shù)分析驗(yàn)證。考慮實(shí)驗(yàn)葉片主要發(fā)生Y方向上的位移，所以提取出實(shí)驗(yàn)葉片在Y方向上的位移曲線。當(dāng)風(fēng)速為8m/s時(shí)，Y方向上各監(jiān)測(cè)點(diǎn)實(shí)驗(yàn)與仿真曲線對(duì)比圖，如圖3所示。

圖3 實(shí)驗(yàn)與仿真曲線對(duì)比圖Fig.3 Comparison of Experimental and Simulation Curves

根據(jù)圖3對(duì)比Y方向上不同監(jiān)測(cè)點(diǎn)曲線圖我們可知：當(dāng)風(fēng)速為8m/s時(shí)，實(shí)驗(yàn)曲線繞仿真曲線呈現(xiàn)上下波動(dòng)，實(shí)驗(yàn)穩(wěn)定值大于仿真穩(wěn)定值。

出現(xiàn)以上問題為葉片實(shí)驗(yàn)環(huán)境激勵(lì)成分復(fù)雜，葉片支架對(duì)監(jiān)測(cè)點(diǎn)位移產(chǎn)生影響等，使實(shí)驗(yàn)極大值和極小值分別大于和小于仿真曲線極值。

當(dāng)風(fēng)速為8m/s時(shí)，Y方向上各監(jiān)測(cè)點(diǎn)實(shí)驗(yàn)與仿真最大值，如表1所示。

表1 實(shí)驗(yàn)與仿真各監(jiān)測(cè)點(diǎn)最大值Tab.1 Maximum Value of Experiment and Simulation at Each Monitoring Point

根據(jù)表1 可知監(jiān)測(cè)點(diǎn)一、二、三實(shí)驗(yàn)與仿真差值分別為0.3098×10-6m、0.4992×10-6m、0.5091×10-5m，所得差值較小，且實(shí)驗(yàn)曲線與仿真曲線波動(dòng)趨勢(shì)一致，可得當(dāng)前仿真模型結(jié)果有效。

5 靜止葉片雙向流固耦合分析

這里使用ANSYS Workbench中Fluent、Transient Structure及System Coupling 三個(gè)模塊進(jìn)行數(shù)值模擬計(jì)算，調(diào)用系統(tǒng)耦合器使數(shù)據(jù)在流體場(chǎng)及固體場(chǎng)之間相互迭代計(jì)算，完成對(duì)葉片的雙向流固耦合計(jì)算，模塊連接方式，如圖4所示。該論文研究的葉片是參考某1000W 垂直軸風(fēng)力機(jī)葉片，三維模型與實(shí)體葉片比例為1：1。

圖4 雙向流固耦合模塊示意圖Fig.4 Schematic Diagram of Two-Way Fluid-Structure Coupling Module

根據(jù)所繪制葉片的三維模型，對(duì)葉片外流場(chǎng)進(jìn)行了針對(duì)性設(shè)計(jì)，對(duì)流場(chǎng)劃分網(wǎng)格后，對(duì)耦合面進(jìn)行優(yōu)化。選中耦合面插入細(xì)節(jié)尺寸Sizing，對(duì)Element size進(jìn)行控制以優(yōu)化耦合面，優(yōu)化后耦合面，如圖5（b）所示。對(duì)固體場(chǎng)外表面網(wǎng)格劃分精度與流體場(chǎng)內(nèi)表面網(wǎng)格劃分精度一致，保證兩側(cè)網(wǎng)格對(duì)應(yīng)關(guān)系更好。由于本論文只針對(duì)垂直軸風(fēng)力機(jī)葉片進(jìn)行流固耦合分析，根據(jù)葉片本身實(shí)際尺寸情況，在制定計(jì)算域時(shí)將其設(shè)置為葉片的數(shù)倍以上。設(shè)置流體范圍為（12×2.4×2.4）m，葉片置于距流場(chǎng)入口3m 處中間位置，葉片弦向與風(fēng)向垂直。所畫計(jì)算域網(wǎng)格數(shù)為1217788，計(jì)算域邊界名稱及其網(wǎng)格模型，如圖5所示。

圖5 計(jì)算域及其網(wǎng)格模型Fig.5 Computational Domain and its Grid Model

5.1 流體場(chǎng)條件分析設(shè)置

對(duì)三維葉片模型處理完成后，在Fluent中對(duì)流場(chǎng)進(jìn)行設(shè)置。該論文對(duì)計(jì)算域湍流模型的選用為：realizable k-ε模型，設(shè)置湍流強(qiáng)度為5%。一般垂直軸風(fēng)機(jī)有效風(fēng)速介于（3～20）m/s之間，適宜風(fēng)速為（6～8）m/s，該論文主要考慮中高風(fēng)速下葉片各位置位移變形情況，所以選取入口風(fēng)速值依次為8m/s、12m/s、15m/s。

計(jì)算前對(duì)流體場(chǎng)各個(gè)條件進(jìn)行設(shè)置，為使本次計(jì)算模擬更加接近風(fēng)機(jī)葉片實(shí)際工作情況，我們認(rèn)為葉片已經(jīng)安裝在風(fēng)機(jī)支架上，即流場(chǎng)四周均不受氣流擾動(dòng)影響，所以流場(chǎng)上下左右四個(gè)面均為對(duì)稱面，前后兩端分別為入口及出口，出口設(shè)置為自由流出，將流場(chǎng)內(nèi)表面定義為雙向流固耦合壁面。雙向流固耦合中固體場(chǎng)會(huì)發(fā)生變化，進(jìn)而影響流場(chǎng)，所以在計(jì)算仿真中采用動(dòng)網(wǎng)格進(jìn)行仿真計(jì)算。最后，合理的時(shí)間步長(zhǎng)對(duì)于求解精度以及對(duì)計(jì)算機(jī)資源的利用可以達(dá)到預(yù)期要求，為求解準(zhǔn)確性以及計(jì)算機(jī)資源的充分使用，設(shè)置時(shí)間步長(zhǎng)為0.005s，每一步保存一次。當(dāng)殘差值達(dá)到1.0×10-3收斂。

5.2 瞬態(tài)結(jié)構(gòu)計(jì)算分析設(shè)置

對(duì)流體場(chǎng)條件設(shè)置完成后，進(jìn)入Transient Structure中對(duì)葉片進(jìn)行設(shè)定。確定葉片材料為6061t6鋁。對(duì)葉片進(jìn)行網(wǎng)格劃分后，將葉片表面做細(xì)化處理，精度與流場(chǎng)內(nèi)表面一致，使結(jié)果更準(zhǔn)確。定義約束面1為遠(yuǎn)端位移，使約束面1沿葉展方向自由移動(dòng)，其余方向全部固定，約束面2為固定約束。在瞬時(shí)分析中的時(shí)間步長(zhǎng)采取與流體場(chǎng)中一致，在此不單獨(dú)設(shè)置，即為0.005s。

6 計(jì)算結(jié)果分析

在對(duì)垂直軸風(fēng)力機(jī)葉片模型進(jìn)行非穩(wěn)態(tài)分析前進(jìn)入System Coupling中對(duì)時(shí)間步長(zhǎng)及截止時(shí)間進(jìn)行合理設(shè)定，保證工作站資源的合理有效利用。考慮到葉片結(jié)構(gòu)為上下對(duì)稱且最大變形出現(xiàn)在靠近尾緣處，我們?cè)谌~片端部尾緣、0.25h尾緣與0.5h尾緣處分別選取一個(gè)監(jiān)測(cè)點(diǎn)，定義為監(jiān)測(cè)點(diǎn)一、二、三，以便觀察葉片在受到流場(chǎng)壓力過程中3個(gè)位置的網(wǎng)格位移變化情況。觀察點(diǎn)位置，如圖1所示。

圖6（a）為當(dāng)來流條件為8m/s時(shí)流體場(chǎng)葉片位置的風(fēng)速流線圖，圖6（b）、圖6（c）為葉片在來流方向及背流方向的壓力分布。由圖6（b）、圖6（c）可知，葉片來流方向壓力大于葉片背流方向壓力。當(dāng)來流條件設(shè)置為12m/s與15m/s時(shí)，壓力變化趨勢(shì)與8m/s時(shí)相同。不同風(fēng)速下葉片位移變化曲線及3個(gè)監(jiān)測(cè)點(diǎn)的總網(wǎng)格位移曲線，如圖7所示。

圖6 風(fēng)場(chǎng)風(fēng)速流線及葉片來流方向、背流方向壓力圖Fig.6 Wind Speed Streamline of Wind Field and Pressure Diagram of Flow Direction and Back Flow Direction of Blade

圖7 不同風(fēng)速下葉片位移變化曲線及3個(gè)監(jiān)測(cè)點(diǎn)的總網(wǎng)格位移曲線Fig.7 The Blade Displacement Curve Under Different Wind Speeds and the Total Grid Displacement Curve of Three Monitoring Points

由圖7 可知：當(dāng)風(fēng)速為8m/s 時(shí)，網(wǎng)格位移量最小為1.5129×10-23m，最大為7.0861×10-6m，當(dāng)風(fēng)速為12m/s 時(shí)，網(wǎng)格位移量最小為1.1328×10-22m，最大為1.1328×10-5m。當(dāng)風(fēng)速為15m/s 時(shí)，網(wǎng)格位移量最小為8.4104×10-19m，最大為2.3977×10-5。

比較3種不同風(fēng)速條件下網(wǎng)格位移量可得，隨著風(fēng)速的增加位移量隨之增加，由監(jiān)測(cè)點(diǎn)一到監(jiān)測(cè)點(diǎn)二到監(jiān)測(cè)點(diǎn)三，最大位移量由大減小再變大，得到葉片產(chǎn)生最大位移的位置。即靠近葉片中部尾緣處。

葉片在受力過程中總網(wǎng)格位移曲線并沒有呈現(xiàn)明顯的上下波動(dòng)，只出現(xiàn)一個(gè)波峰之后便趨于平緩，這是由于仿真葉片在約束選取上并不是使用葉片兩端進(jìn)行約束，而是分別選取了距離中部330mm處與葉片支架連接處進(jìn)行約束。

監(jiān)測(cè)點(diǎn)一、二、三均是在35步左右趨于平緩，可知葉片在受到流場(chǎng)傳遞的壓力時(shí)，中間部位與兩端的振動(dòng)頻率是相等的，網(wǎng)格位移量不同，在靠近葉片中部尾緣部位網(wǎng)格量達(dá)到最大值。

當(dāng)設(shè)置來流條件為8m/s時(shí)監(jiān)測(cè)點(diǎn)一、二、三在0步時(shí)的初始值與達(dá)到穩(wěn)定后的值不同，監(jiān)測(cè)點(diǎn)一的位移增量為2.11×10-6m，監(jiān)測(cè)點(diǎn)二的位移增量為1.30×10-6m，監(jiān)測(cè)點(diǎn)三的位移增量為2.21×10-6m。

0.5h尾緣處網(wǎng)格位移量大于葉片端部尾緣網(wǎng)格位移量大于0.25h尾緣處網(wǎng)格位移量，葉片端部尾緣與0.5h尾緣處網(wǎng)格位移量相差較小，說明葉片主要變形區(qū)域發(fā)生在靠近葉片中部尾緣處。當(dāng)來流條件為12m/s、15m/s時(shí)，位移增量的變化趨勢(shì)與來流條件8m/s時(shí)一致。

不同風(fēng)速下各個(gè)監(jiān)測(cè)點(diǎn)的曲線對(duì)比圖，如圖8所示。由圖8可知，同一監(jiān)測(cè)點(diǎn)在不同風(fēng)速下初始值也不同，隨著風(fēng)速的增加，總位移初始值隨之增加，穩(wěn)定值也越大。

圖8 不同風(fēng)速下各個(gè)監(jiān)測(cè)點(diǎn)曲線對(duì)比圖Fig.8 Curve Comparison of Each Monitoring Point Under Different Wind Speeds

經(jīng)過計(jì)算葉片端部尾緣的位移增量，可知不同來流情況下該點(diǎn)的位移增量值：風(fēng)速為8m/s時(shí)，位移增量為2.11×10-6m；風(fēng)速為12m/s 時(shí)，位移增量為3.36×10-6m；風(fēng)速為15m/s 時(shí)，位移增量為4.88×10-6m。

當(dāng)來流條件為15m/s 時(shí)，葉片端部尾緣處總網(wǎng)格位移增量最大。

由圖可知，隨著風(fēng)速的增加，葉片端部網(wǎng)格位移增加量隨之增加，監(jiān)測(cè)點(diǎn)曲線波動(dòng)幅值隨之增加。0.25h尾緣與0.5h尾緣處監(jiān)測(cè)點(diǎn)對(duì)比曲線有同樣的趨勢(shì)。

沿葉片展向各方向與總位移的曲線對(duì)比圖，如圖9所示。考慮到葉片上下對(duì)稱，葉片取值范圍為（0～0.5）h。

圖9 沿葉片展向各方向與總位移的曲線對(duì)比圖Fig.9 Comparison Diagram of Displacement Curves Along Each Direction of Blade Spanning

由圖9可知葉片主要變形出現(xiàn)在（0.25～0.35）h處，沿各個(gè)方向上的位移最大值均出現(xiàn)在（0.35～0.4）h處。

沿著葉片展向x，y，z三個(gè)方向上的位移先增加后減少，然后再增加減少，另一半葉片曲線變化趨勢(shì)與此相同。

y方向上的位移大于x方向與z方向的位移，并與總位移的趨勢(shì)相同。

對(duì)比0.5h處的位移值，當(dāng)風(fēng)速為8m/s時(shí)，y方向的位移值為3.44×10-6m，而總位移的位移值為3.47×10-6m，兩者相差不大，即葉片在受力過程中主要發(fā)生y方向上的變形。

對(duì)比葉片y方向上的位移量，沿葉片展向位移量呈現(xiàn)上下波動(dòng)，伴隨風(fēng)速的增加位移量也增加，可知風(fēng)速與葉片在耦合面上的擺振成正比。

7 結(jié)論

該論文合理應(yīng)用ANSYS Workbench實(shí)現(xiàn)垂直軸風(fēng)力機(jī)系統(tǒng)中葉片的雙向流固耦合作用，使用實(shí)驗(yàn)和仿真相互驗(yàn)證的方法，對(duì)葉片進(jìn)行合理分析。根據(jù)計(jì)算結(jié)果對(duì)葉片的位移變形進(jìn)行分析，分析結(jié)果表明：（1）隨著風(fēng)速的增加，葉片位移增加明顯，0.5h處大于葉片端點(diǎn)處大于0.25h處。（2）風(fēng)載與葉片位移成正比，最大變形出現(xiàn)在（0.25～0.35）h、（0.65～0.75）h 處，最大位移出現(xiàn)在（0.35～0.4）h、（0.6～0.65）h處，y軸正方向與來流方向相同，風(fēng)載是影響葉片產(chǎn)生變形的主要載荷，且沿葉片展向位移值上下波動(dòng)。（3）通過開展風(fēng)力機(jī)葉片流場(chǎng)和結(jié)構(gòu)場(chǎng)雙向流固耦合數(shù)值模擬計(jì)算可知，風(fēng)載與葉片在耦合面上的擺振成正比，且沿葉片展向上下波動(dòng)。