汪 泉，陳曉田，黃 攀，甘 笛，曾利磊，胡夢(mèng)杰

（湖北工業(yè)大學(xué) 機(jī)械工程學(xué)院，湖北 武漢 430068）

0 前言

風(fēng)能作為一種無(wú)污染、 無(wú)危害并且儲(chǔ)量十分豐富的可再生能源，越來(lái)越受到各國(guó)學(xué)者的青睞。與水平軸風(fēng)力機(jī)相比，H 型垂直軸風(fēng)力機(jī)結(jié)構(gòu)簡(jiǎn)單，不須要設(shè)置偏航裝置，發(fā)電機(jī)安裝在塔底，便于安裝和維修，而且葉片的制造過(guò)程簡(jiǎn)單，不需要較高的成本[1]。

近年來(lái)研究人員對(duì)H 型垂直軸風(fēng)力機(jī)的葉片進(jìn)行了針對(duì)性研究。 Ning Ma[2]基于遺傳算法對(duì)NACA0018 翼型進(jìn)行優(yōu)化，CFD 仿真結(jié)果顯示：當(dāng)葉尖速比（λ）為0.9 時(shí)，較優(yōu)化前風(fēng)力機(jī)功率系數(shù)（Cp）提高了 26.82%。Changping Liang[3]采用類函數(shù)和形狀函數(shù)變參數(shù)優(yōu)化方法對(duì)NACA0015 翼型幾何形狀進(jìn)行擾動(dòng)，最大升力系數(shù)、最大升阻比和最大切向力系數(shù) （CT） 分別提高 7.5%，9%和8.87%。 梁昌平[4]利用CST 參數(shù)化方法對(duì)NACA0015翼型幾何外形進(jìn)行擾動(dòng)，基于非支配排序的遺傳算法進(jìn)行幾何優(yōu)化，優(yōu)化后的翼型最大升力系數(shù)提高 6.63%，最大 CT 提高 22.7%。 Carlos Sim?o Ferreira[5]推導(dǎo)垂直軸風(fēng)力機(jī)翼型氣動(dòng)優(yōu)化公式，并驗(yàn)證了優(yōu)化公式的有效性。

上述研究主要以對(duì)稱翼型為研究對(duì)象，而對(duì)非對(duì)稱翼型的垂直軸風(fēng)力機(jī)研究頗少[6]。 葉片作為風(fēng)力機(jī)最重要的組成部分，是捕獲風(fēng)能的重要部件[7]。 在不同的風(fēng)能資源情況下，翼型的選擇對(duì)提高H 型垂直軸風(fēng)力機(jī)的風(fēng)能利用率有著重要的作用。 本文基于 NACA0018，NACA4418 翼型，以風(fēng)力機(jī)的Cp和單個(gè)葉片的CT 作為判斷依據(jù)，建立垂直軸風(fēng)力機(jī)模型并進(jìn)行CFD 仿真，對(duì)比分析在不同的λ 下兩種風(fēng)力機(jī)模型的氣動(dòng)特性。

1 CFD模型建立

1.1 幾何模型與邊界條件設(shè)置

NACA0018 和 NACA4418 翼型如圖 1 所示。

圖 1 NACA0018 和 NACA4418 翼型Fig.1 NACA0018 airfoil and NACA4418 airfoil

本文中兩種翼型的弦長(zhǎng)均為0.1 m。兩種垂直軸風(fēng)力機(jī)模型的主要參數(shù)如表1 所示。

表1 風(fēng)力機(jī)模型參數(shù)Table 1 Model parameters of wind turbine

連桿作為風(fēng)力機(jī)的一部分，同樣受到空氣阻力。當(dāng)連桿的直徑很小時(shí)，對(duì)風(fēng)力機(jī)的功率轉(zhuǎn)換影響微乎其微，故模型中省略連桿。

圖2 所示的簡(jiǎn)化幾何模型與實(shí)驗(yàn)?zāi)Ｐ拖嗨芠8]。矩形的DE 邊為速度入口，入口邊界條件設(shè)為速度入口，速度大小為4.1 m/s，方向從左到右。 出口邊界為矩形的FG 邊，出口邊界條件設(shè)為壓力出口，為標(biāo)準(zhǔn)大氣壓。矩形計(jì)算域的上下邊界條件設(shè)為對(duì)稱。 圖2 中的 C 區(qū)域?yàn)殪o止區(qū)域，A 和 B 區(qū)域?yàn)樾D(zhuǎn)區(qū)域，不同區(qū)域之間的交界面設(shè)為interface。 旋轉(zhuǎn)區(qū)域的滑移網(wǎng)格運(yùn)動(dòng)方式設(shè)為moving mesh。 翼型壁面設(shè)為無(wú)滑移壁面條件。

圖2 CFD 計(jì)算模型Fig.2 Computational model of CFD

1.2 網(wǎng)格劃分與計(jì)算模型時(shí)間步長(zhǎng)設(shè)置

對(duì)于CFD 分析，一般網(wǎng)格劃分占整個(gè)前處理周期的80%，網(wǎng)格的質(zhì)量對(duì)計(jì)算結(jié)果和計(jì)算精度有很大的影響。 結(jié)構(gòu)化網(wǎng)格的質(zhì)量要高于非結(jié)構(gòu)網(wǎng)格的質(zhì)量，計(jì)算精度和計(jì)算效率都比非結(jié)構(gòu)化網(wǎng)格好。為保證仿真結(jié)果的準(zhǔn)確性，本文網(wǎng)格劃分均采用結(jié)構(gòu)化網(wǎng)格，采用CFD 前處理工具ICEM進(jìn)行網(wǎng)格劃分，網(wǎng)格劃分如圖3 所示。

圖3 網(wǎng)格劃分Fig.3 Mesh generation

以NACA0018 翼型建立的二維模型網(wǎng)格數(shù)量為470 258，節(jié)點(diǎn)數(shù)為 466 585。 以翼型 NACA4418建立的二維模型網(wǎng)格數(shù)量為563 138，節(jié)點(diǎn)數(shù)為559 465。 由于近翼型壁面處的流場(chǎng)復(fù)雜，考慮葉片近壁面處網(wǎng)格質(zhì)量對(duì)計(jì)算結(jié)果有影響，臨近翼型壁面處的網(wǎng)格均加密處理。

為了保證計(jì)算結(jié)果的精度和可信度，本文對(duì)H 型垂直軸風(fēng)力機(jī)的研究采用Fluent 滑移網(wǎng)格技術(shù)。 時(shí)間步長(zhǎng)的選擇關(guān)系到計(jì)算的效率和結(jié)果的準(zhǔn)確性。 在Fluent 軟件中預(yù)先設(shè)置計(jì)算的時(shí)間步數(shù)和時(shí)間步長(zhǎng)。 時(shí)間步長(zhǎng)越小，計(jì)算精度就越高，反之計(jì)算精度就越低。 在設(shè)定好一個(gè)時(shí)間步長(zhǎng)內(nèi)風(fēng)力機(jī)旋轉(zhuǎn)的度數(shù)后，得到風(fēng)力機(jī)在旋轉(zhuǎn)一個(gè)周期內(nèi)的步數(shù)（K）和旋轉(zhuǎn) n 度的時(shí)間步長(zhǎng)（t）。

式中：ω 為風(fēng)力機(jī)轉(zhuǎn)速。

1.3 湍流模型

Fluent 提供的湍流模型很多，每一種湍流模型的適用情形都各不相同。 單方程Spalart-Allmaras 適用于不是特別復(fù)雜壓力梯度作用下的邊界層流動(dòng)。兩方程主要有k-epslion 和k-omega 兩種湍流模型。標(biāo)準(zhǔn)的k-epslion 湍流模型具有較好的通用性，但遇到旋渦以及逆壓梯度問(wèn)題時(shí)，計(jì)算精度較低。 k-omega 湍流模型能較好地解決邊界層、自由剪切的問(wèn)題。特別是SST k-omega 模型適用于逆壓梯度下復(fù)雜邊界層流動(dòng)及分離。 梅毅[9]分別用RNG k-epslion 模型、Realizable k-epslion模型和SST k-omega 模型進(jìn)行垂直軸風(fēng)力機(jī)的數(shù)值模擬，結(jié)果表明SST k-omega 模型計(jì)算得到的功率系數(shù)更接近實(shí)驗(yàn)值。綜上所述，基于本文中建立的垂直軸風(fēng)力機(jī)二維流場(chǎng)模型為瞬態(tài)，且風(fēng)輪附近流場(chǎng)較復(fù)雜，故選擇SST k-omega 湍流模型。

1.4 求解器設(shè)置

Fluent 軟件中提供了壓力基求解器和密度基求解器兩種求解方式。 本文選擇壓力基求解器。SIMPLE 算法和SIMPLEC 算法都適用于大多數(shù)不可壓縮流動(dòng)的工程問(wèn)題，但SIMPLEC 算法在四邊形網(wǎng)格的計(jì)算中更準(zhǔn)確，收斂速度更快，故采用SIMPLEC 算法。 壓力、動(dòng)量和其他變量分別采用二階逆風(fēng)、有界中心差分和三階MUSCL 格式，各項(xiàng)殘差為0.000 01。

2 計(jì)算結(jié)果處理方法

在Fluent 仿真模擬中，監(jiān)測(cè)到的圖像為風(fēng)力機(jī)旋轉(zhuǎn)過(guò)程中的力矩系數(shù)曲線，力矩系數(shù)輸出值的處理很重要。 在仿真模擬的前幾個(gè)周期由于風(fēng)力機(jī)處于非穩(wěn)定狀態(tài)，監(jiān)測(cè)到的力矩系數(shù)曲線雜亂無(wú)章[9]。 待風(fēng)力機(jī)穩(wěn)定后，力矩系數(shù)曲線呈周期性變化。為保證計(jì)算結(jié)果的準(zhǔn)確性，待風(fēng)力機(jī)處于穩(wěn)定狀態(tài)過(guò)后，取最后一個(gè)周期的力矩系數(shù)數(shù)據(jù)進(jìn)行分析。 由Fluent 仿真模擬中輸出的力矩系數(shù)求得扭矩（T）。

式中：ρ 為空氣密度，g/m3；V∞為來(lái)流風(fēng)速，m/s；L為葉片弦長(zhǎng)，m；Cm為力矩系數(shù)。

由扭矩可以求得風(fēng)能利用系數(shù)（Cp）。

式中：A 為風(fēng)力機(jī)掃風(fēng)面積，m2。

3 結(jié)果與分析

3.1 仿真結(jié)果驗(yàn)證

以模型1 為例，與風(fēng)洞實(shí)驗(yàn)中建立的模型相似。利用風(fēng)洞實(shí)驗(yàn)結(jié)果與CFD 仿真模擬結(jié)果作對(duì)比，驗(yàn)證CFD 仿真結(jié)果的正確性。 圖4 為風(fēng)洞實(shí)驗(yàn)與CFD 仿真在不同λ 下垂直軸風(fēng)力機(jī)的Cp。由圖4 可知，CFD 仿真結(jié)果較好地顯示了Cp的變化趨勢(shì)。 在低λ 范圍內(nèi)，風(fēng)洞實(shí)驗(yàn)與CFD 仿真結(jié)果基本一致。 在高λ 范圍內(nèi)，CFD 仿真模擬低估了H 型垂直軸風(fēng)力機(jī)的Cp，但仿真結(jié)果在可接受的誤差范圍內(nèi)。 從Cp整體的變化趨勢(shì)來(lái)看，CFD仿真模擬與風(fēng)洞實(shí)驗(yàn)結(jié)果有較好的吻合性。

圖4 風(fēng)洞實(shí)驗(yàn)與CFD 仿真模擬在不同λ 下的功率系數(shù)Fig.4 Power coefficient of numerical simulation and wind tunnel test at different tip seed ratio

3.2 風(fēng)能利用率

圖5 不同λ 下的CpFig.5 Power coefficient at different tip seed ratio

取風(fēng)力機(jī)仿真結(jié)果中最后一個(gè)周期的力矩系數(shù)的平均值，根據(jù)式（1），（2）求出風(fēng)力機(jī)的 Cp。 圖5 為兩種模型在不同 λ 下的 Cp。由圖5 可知，兩種模型在不同λ 下的Cp曲線變化趨勢(shì)大致相同。當(dāng)風(fēng)力機(jī)在低速狀態(tài)（低λ）下運(yùn)行時(shí)，風(fēng)力機(jī)Cp較低。 因?yàn)轱L(fēng)力機(jī)在低速狀態(tài)下運(yùn)行時(shí)，葉片多數(shù)處于失速狀態(tài)，產(chǎn)生的氣動(dòng)力較小，因此，風(fēng)力機(jī)氣動(dòng)特性較差。 當(dāng)風(fēng)力機(jī)在高速狀態(tài)（高λ）下運(yùn)行時(shí)，風(fēng)力機(jī)Cp隨轉(zhuǎn)速增加而降低，因?yàn)楦咚贍顟B(tài)時(shí)，風(fēng)輪局部氣流受到破壞，導(dǎo)致風(fēng)力機(jī)氣動(dòng)特性差。 當(dāng)風(fēng)力機(jī)在最佳λ 狀態(tài)運(yùn)行時(shí)，風(fēng)力機(jī)具有良好的氣動(dòng)性能，Cp達(dá)到最大值。 此外，在低 λ 時(shí)，NACA448 翼型的風(fēng)力機(jī) Cp較 NACA0018翼型的風(fēng)力機(jī)高；相反，在高 λ 時(shí)，NACA0018 翼型的風(fēng)力機(jī)Cp高。說(shuō)明NACA4418 翼型風(fēng)力機(jī)在低速狀態(tài)下的氣動(dòng)特性好，適用于低風(fēng)速地區(qū)。而NACA0018 翼型風(fēng)力機(jī)在高速狀態(tài)下的氣動(dòng)特性好，適用于高風(fēng)速地區(qū)。

3.3 不同λ下翼型的CT分析

取風(fēng)力機(jī)最后一個(gè)周期的單個(gè)葉片的切向力計(jì)算結(jié)果作為分析對(duì)象。 圖6 為單個(gè)翼型在低λ下的 CT 變化曲線。 由圖6 可知： 在低的 λ 范圍內(nèi)，由于整個(gè)風(fēng)力機(jī)處于失速狀態(tài)，葉片受到失速氣流的干擾，導(dǎo)致CT 曲線波動(dòng)雜亂無(wú)章；在整個(gè)方位角范圍內(nèi)，NACA0018 翼型CT 曲線波動(dòng)范圍比NACA4418 翼型大；在方位角180°附近，當(dāng)λ為 1.53，1.93 和 2.46 時(shí)，NACA4418 翼型 CT 大得多；當(dāng) λ 為 1.53 和 1.93 時(shí)，在方位角 300 °附近，NACA0018 翼型的負(fù)CT 值同樣比NACA0018 翼型大得多。 這也是在低λ 時(shí)NACA4418 翼型風(fēng)力機(jī)風(fēng)能利用率較NACA0018 翼型的高的原因。 由此可知，相比較 NACA0018 翼型，NACA4418 翼型在低風(fēng)速下有更好的氣動(dòng)性能。

圖6 單個(gè)葉片在低λ 下的CT 曲線Fig.6 Tangential force coefficient curves of single blade at low tip speed ratio

圖7 為單個(gè) NACA0018 與 NACA4418 翼型在高λ 下CT 變化曲線。

圖7 高λ 下單個(gè)葉片的CT 曲線Fig.7 Tangential force coefficient curves of single blade at low tip speed ratio

由圖7 可知，NACA0018 與 NACA4418 翼型CT 曲線變化趨勢(shì)相同。隨著λ 的增大，翼型CT曲線下移，風(fēng)力機(jī)Cp隨著轉(zhuǎn)速的增加逐漸降低，說(shuō)明翼型在高λ 下，葉片的氣動(dòng)力下降。 翼型方位角在 0～180 °時(shí)，NACA0018 與 NACA4418 翼型CT曲線的波峰值隨著轉(zhuǎn)速的增加差值增大，且NACA0018 翼型的CT 平均值高于 NACA4418 翼型，說(shuō)明NACA0018 翼型在上風(fēng)區(qū)氣動(dòng)性能較NACA4418 翼型好。 而翼型在方位角為180～360°時(shí)，翼型CT 隨著λ 的增加向縱軸負(fù)方向移動(dòng)，但NACA4418 翼型 CT 曲線大部分在 NACA0018 翼型的上方，且NACA4418 翼型的CT 平均值高于NACA0018 翼型，說(shuō)明 NACA4418 翼型在下風(fēng)區(qū)的氣動(dòng)性能較NACA0018 翼型好。

4 結(jié)論

本文通過(guò)對(duì)垂直軸風(fēng)力機(jī)兩種翼型氣動(dòng)特性的比較研究，得到以下結(jié)論。

①當(dāng)垂直軸風(fēng)力機(jī)翼型為NACA4418 翼型，在低λ 的情況下運(yùn)行時(shí)，風(fēng)力機(jī)的氣動(dòng)特性要優(yōu)于 NACA0018 翼型的風(fēng)力機(jī)。 相反，在高的 λ 時(shí)，NACA0018 翼型風(fēng)力機(jī)的氣動(dòng)特性好。 因此，NACA4418 翼型適用于低風(fēng)速區(qū)域，而NACA0018 翼型適用于高風(fēng)速區(qū)域。

②風(fēng)力機(jī)在高λ 時(shí)：NACA0018 翼型在上風(fēng)區(qū)的CT 平均值高于NACA4418 翼型；在下風(fēng)區(qū)，NACA4418 翼型的 CT 平均值高于 NACA0018 翼型。說(shuō)明高λ 時(shí)，NACA0018 翼型在上風(fēng)區(qū)的氣動(dòng)性能比NACA4418 翼型好，在下風(fēng)區(qū)NACA4418翼型氣動(dòng)性能好。