， ， ，，

(1.中海油研究總院，北京 100027;

2.中國石油大學(xué)博士后流動站，北京 102200;

3.中國石油大學(xué)機(jī)械與儲運(yùn)工程學(xué)院，北京 102200;

4. 哈爾濱工程大學(xué)海洋可再生能源研究所，哈爾濱 150001)

0 引 言

潮流能作為一種新興的清潔、可再生能源正在為世界各種爭相開發(fā)，而垂直軸水輪機(jī)作為潮流能轉(zhuǎn)換的核心裝置之一，已成為潮流能開發(fā)的主要途徑。水輪機(jī)技術(shù)的研究主要來自實(shí)驗(yàn)驗(yàn)證與計算機(jī)的數(shù)值模擬。垂直軸水輪機(jī)和垂直軸風(fēng)力機(jī)基本原理是相同的，在數(shù)值模擬方面可以相互借鑒。2005年，Horiuchi等人[1]使用STAR-CD 模擬二維假設(shè)下的垂直軸風(fēng)力機(jī)速度場，并與試驗(yàn)值進(jìn)行比較。2007年，Guerri等人[2]使用STAR-CD模擬二維垂直軸風(fēng)力機(jī)葉片的受力。2008年孫科提出應(yīng)用滑移網(wǎng)格理論研究H型葉輪的動邊界CFD數(shù)值模擬問題[3-4]。2009年，Hwang等人[5]應(yīng)用STAR-CD，采用k-ε湍流模型對擺線式垂直水輪機(jī)進(jìn)行了研究。

CFD數(shù)值模擬的方法在水輪機(jī)技術(shù)研究過程中得到了廣泛的應(yīng)用，但是模擬結(jié)果的精度及可信度一直是所有研究者共同關(guān)注的焦點(diǎn)。盡管垂直軸水輪機(jī)或風(fēng)力機(jī)CFD數(shù)值預(yù)報方面取得了很大的發(fā)展，但是有限的精度仍舊迫切需要得到改進(jìn)。這其中如物理模型、網(wǎng)格劃分、邊界條件、湍流模型、CFD求解器等因素引起的計算結(jié)果中的誤差，都是不能被忽略的。因此，本文針對CFD方法在垂直軸水輪機(jī)數(shù)值模擬中的有效性進(jìn)行了詳細(xì)的研究，包括湍流模型、網(wǎng)格尺度，以及時間步長對計算精度的影響。

1 CFD理論基礎(chǔ)

計算流體力學(xué)(Computational Fluid Dynamics，簡稱CFD)是流體力學(xué)研究領(lǐng)域中的重要技術(shù)之一，其基本思路為使用數(shù)值方法在計算機(jī)中對流體力學(xué)的控制方程進(jìn)行求解，從而達(dá)到模擬預(yù)測流場的目的。下面介紹本文中涉及的湍流模型及其特點(diǎn)，公式的推導(dǎo)過程可以參考Ansys CFX用戶使用手冊[6]。

1.1 k-ε模型

(1)標(biāo)準(zhǔn)的k-ε模型：

標(biāo)準(zhǔn)k-ε模型：是從實(shí)驗(yàn)現(xiàn)象總結(jié)出來的半經(jīng)驗(yàn)公式，也是最簡單的兩方程模型，主要基于湍流動能k和湍流耗散率ε。k方程是個精確方程，ε方程是個由經(jīng)驗(yàn)公式導(dǎo)出的方程。該模型假定流場完全是湍流，分子之間的粘性可以忽略。因此，只對完全是湍流的流場有效。

(2)RNGk-ε模型：

1.2 k-ω模型

(1)標(biāo)準(zhǔn)的k-ω模型：

兩方程模型，包含了低雷諾數(shù)影響、可壓縮性影響和剪切流擴(kuò)散，因此適用于尾跡流動計算、混合層計算、射流計算以及受壁面限制的流動計算和自由剪切流計算。

(2)SSTk-ω模型：

它綜合了k-ω模型在近壁區(qū)計算的優(yōu)點(diǎn)和k-ε模型在遠(yuǎn)場計算的優(yōu)點(diǎn)，將k-ω模型和標(biāo)準(zhǔn)k-ε模型都乘以一個混合函數(shù)后再相加就得到這個模型。與標(biāo)準(zhǔn)k-ω相比，SSTk-ω模型中增加了橫向耗散導(dǎo)數(shù)項，同時在湍流粘性定義中考慮了湍流剪切應(yīng)力的輸運(yùn)過程，模型中使用的湍流常數(shù)也有所不同。這些特點(diǎn)使得 SSTk-ω模型的適用范圍更廣，比如可以用于帶逆壓梯度的流動計算、翼型計算、跨音速激波流動計算等等。

1.3 雷諾應(yīng)力模型(RSM)

雷諾應(yīng)力模型中沒有渦粘度的各向同性假設(shè)，因此從理論上來說湍流模式理論要更精確。從理論上說，RSM 應(yīng)該比一方程和兩方程模型的計算精度更高，但實(shí)際上雷諾應(yīng)力模型的精度受限于模型的封閉形式，因此雷諾應(yīng)力明顯具有各向異性的特點(diǎn)時才必須使用 RSM，比如龍卷風(fēng)、燃燒室內(nèi)流動等帶有強(qiáng)烈旋轉(zhuǎn)的流動問題。

2 垂直軸水輪機(jī)的CFD有效性驗(yàn)證

2.1 計算算例

1983年Strikland等人[5]等人做的一系列研究垂直軸風(fēng)力機(jī)動態(tài)失速試驗(yàn)中，測量和記錄的數(shù)據(jù)包括葉片的瞬時切向和法向受力、葉片表面壓力分布以及尾流流場等，試驗(yàn)結(jié)果作為一個可靠的參考依據(jù)。 Deglaire等人[7]、Hwang等人、Tchon和Paraschivoiu[8]和汪魯兵[9]的模擬研究中都是與其比較得出結(jié)論的。本文研究同樣采用Strikland等人的試驗(yàn)?zāi)Ｐ蜑轵?yàn)證算例。模型的具體參數(shù)如表1所示。

表1 Strickland試驗(yàn)輪機(jī)的主要參數(shù)

2.2 網(wǎng)格劃分

網(wǎng)格劃分是CFD模擬過程中較為耗時的環(huán)節(jié)，也是直接影響模擬精度和效率的關(guān)鍵因素之一。對于所有的CFD數(shù)值模擬，結(jié)構(gòu)網(wǎng)格模擬精度和速度要明顯優(yōu)于非結(jié)構(gòu)網(wǎng)格情況，因此本文系采用結(jié)構(gòu)化網(wǎng)格對流場進(jìn)行劃分。根據(jù)垂直軸水輪機(jī)的結(jié)構(gòu)特點(diǎn)，葉片在展長方向上具有相同的橫截面，為了節(jié)省計算時間，本文簡化為二維模型。

水蛭交配一般發(fā)生在早上或夜間，交配后一個月左右，雌體生殖器分泌出稀薄的黏液，中包被卵帶，形如“蠶繭”，排出體外，在濕泥中孵化。產(chǎn)繭一般在土壤（含30%水分）中進(jìn)行，產(chǎn)繭時間一般為4月中下旬溫度在20-25℃以上。

劃分網(wǎng)格之前，還需要對模型中的流場進(jìn)行基本的預(yù)估。過疏或過密的網(wǎng)格都會對計算有極大影響。過疏的網(wǎng)格往往會得到不精確甚至完全錯誤的結(jié)果；過密的網(wǎng)格會使計算量增大，使計算難以收斂。文中的計算網(wǎng)格如圖1和圖2所示。

圖1 Turbine 網(wǎng)格分布圖

圖2 Turbine 葉片周圍網(wǎng)格分布圖

2.3 邊界條件的設(shè)置

流體域內(nèi)是以大氣壓為參考?xì)鈮旱?；在進(jìn)口邊界處設(shè)置為速度進(jìn)口條件，給定均勻來流速度V及湍流參數(shù)；出口邊界采用開敞邊界；流體計算域的側(cè)面采用滑動墻面(slip wall)設(shè)定，這種墻面設(shè)定可以有效地提高計算速度，葉片表面為不可滑移壁面；旋轉(zhuǎn)域給定旋轉(zhuǎn)角速度，旋轉(zhuǎn)域和靜止域通過交界面連接；使用入口處流動參數(shù)來初始化流場；采用瞬態(tài)求解器。

圖3 邊界條件設(shè)定

3 結(jié)果對比

3.1 不同湍流模型計算結(jié)果對比

圖4給出了標(biāo)準(zhǔn)、RNG 、SST和Omega雷諾應(yīng)力四種湍流模型計算得到的水輪機(jī)在速比為5.0時，葉片法向力系數(shù)和切向力系數(shù)隨位置角的變化歷程與試驗(yàn)測量結(jié)果的比較。從總體來看，不同湍流模型計算結(jié)果和試驗(yàn)結(jié)果趨勢基本一致，法向力系數(shù)要比切向力系數(shù)吻合的好。但是在數(shù)值上，不同湍流模型結(jié)果還是存在一些差異。

從圖4(a)法向力系數(shù)曲線可以看出，結(jié)果的主要差別出現(xiàn)在谷值和峰值上，四種湍流模型模擬結(jié)果谷值均大于試驗(yàn)結(jié)果，雷諾應(yīng)力模型結(jié)果偏差更大一些，其他三種模型結(jié)果相當(dāng)；模擬結(jié)果峰值均小于試驗(yàn)值，SST模型模擬結(jié)果與試驗(yàn)值更為接近，其次是雷諾應(yīng)力模型，標(biāo)準(zhǔn)模型和RNG 模型結(jié)果相當(dāng)。從4(b)切向力系數(shù)曲線可以看出，四種湍流模型模擬結(jié)果在上游盤面()與試驗(yàn)結(jié)果符合程度好于下游盤面時的情況，這主要是因?yàn)橄掠伪P面經(jīng)過葉片擾動，流場變得更加復(fù)雜，使模擬精度下降。在上游盤面谷值處四種湍流模型模擬結(jié)果很接近，均略大于試驗(yàn)結(jié)果，同時與試驗(yàn)結(jié)果相位也有偏差；上游盤面峰值處，RNG 模型結(jié)果高于試驗(yàn)結(jié)果，標(biāo)準(zhǔn)模型結(jié)果低于試驗(yàn)結(jié)果，SST模型和雷諾應(yīng)力模型結(jié)果與試驗(yàn)結(jié)果都比較接近；下游盤面谷值處四種湍流模型模擬結(jié)果很接近，但是都小于試驗(yàn)值，同樣也存在相位偏差；下游盤面峰值情況，四種湍流模型模擬結(jié)果均出現(xiàn)了平峰，這一點(diǎn)與試驗(yàn)值差別較大，SST模型模型結(jié)果更接近試驗(yàn)值。

總的來說，采用SST模型所計算的結(jié)果與試驗(yàn)結(jié)果吻合更好，這也進(jìn)一步體現(xiàn)了SST模型對逆向壓力梯度與流動分離的流場模擬的優(yōu)勢，雖然采用雷諾應(yīng)力模型計算結(jié)果與SST模型計算的準(zhǔn)確性相差不大，但是雷諾應(yīng)力模型比SST模型計算耗時更長，因此在本文以后的計算中選用SST模型。

3.2 網(wǎng)格尺度對計算結(jié)果的影響

網(wǎng)格的疏密程度與計算迭代的收斂速度以及最終的計算精度緊密相關(guān)。網(wǎng)格節(jié)點(diǎn)數(shù)量大則意味著計算結(jié)果精度較高，但這是以增加計算時間為代價的。通常情況下，網(wǎng)格的生成過程中人的主觀因素對網(wǎng)格的質(zhì)量將產(chǎn)生很大的影響。 CFD 用戶在生成網(wǎng)格時僅僅憑個人經(jīng)驗(yàn)，在平衡計算時間和精度的前提下確定網(wǎng)格分布的拓?fù)浣Y(jié)構(gòu)、密度等。文中需要對垂直軸水輪機(jī)非定常繞流流場做大量計算，希望網(wǎng)格在保證一定精度的前提下，網(wǎng)格數(shù)盡可能少。因此在開展數(shù)值計算前，有必要對網(wǎng)格進(jìn)行測試，測試具有不同疏密程度的網(wǎng)格系統(tǒng)對垂直軸水輪機(jī)計算結(jié)果的影響。

垂直軸水輪機(jī)數(shù)值模擬中，實(shí)質(zhì)是葉片的繞流問題，因此葉片周圍網(wǎng)格最為重要，為了進(jìn)行網(wǎng)格測試，葉片周圍生成4套密度不同的“⊙”形結(jié)構(gòu)網(wǎng)格，周向節(jié)點(diǎn)布置相同，僅變化徑向的網(wǎng)格密度，相關(guān)信息見表2。

表2 不同網(wǎng)格結(jié)構(gòu)相關(guān)信息

仍然以水輪機(jī)在速比5.0時的試驗(yàn)結(jié)果與數(shù)值結(jié)果進(jìn)行對比，以此判斷所選網(wǎng)格密度是否合理。圖5給出了四種不同網(wǎng)格密度(Mesh 1～4)計算得到葉片法向力系數(shù)和切向力系數(shù)隨位置角的變化歷程與試驗(yàn)測量結(jié)果的比較。從圖5(a)法向力系數(shù)的對比曲線看出，網(wǎng)格疏密對法向力的影響不大，四種不同網(wǎng)格密度的計算結(jié)果十分接近。而網(wǎng)格疏密對切向力的影響比較大，從圖5(b)切向力系數(shù)的對比曲線看出，隨著網(wǎng)格密度的增加而呈越來越接近試驗(yàn)結(jié)果的趨勢，但是Mesh2和Mesh1的結(jié)果差別不大，基本滿足網(wǎng)格無關(guān)性的要求。

圖5 葉片法向力和切向力系數(shù)隨位置角的變化(λ=5.0)

這也和文獻(xiàn)推薦值很吻合，當(dāng)時可以充分利用低雷諾數(shù)模型預(yù)測分離流的優(yōu)勢。鑒于速比5.0是一個比較高的速比，分離流動還不是很明顯，因此文中后面計算全部采用網(wǎng)格密度類型Mesh2。

3.3 時間步長對計算結(jié)果的影響

在非穩(wěn)態(tài)問題求解時，時間步長影響模擬的收斂性及CPU計算時間。正如存在網(wǎng)格無關(guān)解一樣，非穩(wěn)態(tài)數(shù)值模擬應(yīng)當(dāng)存在時間步長無關(guān)解：數(shù)值解在工程允許的偏差范圍內(nèi)不再隨時間步長的減小而發(fā)生變化。

但是，時間步長的減小導(dǎo)致了計算時間的迅速增加。由于描述計算對象的微分方程的剛性較強(qiáng)，較小的時間步長有助于得到較為準(zhǔn)確的結(jié)果，且能提高時間步內(nèi)的收斂速度；但是若計算的總時間一定，這并不足以抵消時間步數(shù)增加而引起的計算量增加。時間步長往往需要在計算精度和收斂速度之間進(jìn)行折衷。下面采用簡單的增加時間步長的方法來分析時間步長對垂直軸水輪機(jī)數(shù)值模擬結(jié)果的影響，時間步長主要根據(jù)葉輪每步轉(zhuǎn)過角度來設(shè)定，采用上節(jié)中介紹的Mesh2網(wǎng)格模型，時間步長設(shè)置見表3。

表3 Strickland實(shí)驗(yàn)輪機(jī)的主要參數(shù)

圖6給出了四個不同時間步長計算得到的水輪機(jī)在速比為5.0時，葉片法向力系數(shù)和切向力系數(shù)隨位置角的變化歷程與試驗(yàn)測量結(jié)果的比較。從圖中可以看出，時間步長對計算結(jié)果的影響和網(wǎng)格影響很相似，都是對法向力影響比較小，對切向力影響比較大，湍流模型的影響也有類似結(jié)論，究其原因：一是因?yàn)槿~片切向力系數(shù)的大小比法向力系數(shù)小一個量級；二是切向力系數(shù)的變化規(guī)律比法向力系數(shù)復(fù)雜，在一個周期內(nèi)，切向力系數(shù)在上游盤面和下游盤面各出現(xiàn)一次峰值和谷值，即其變化頻率比法向力系數(shù)大。從圖6(b)葉片切向力系數(shù)曲線中可以看出，時間步長越大，計算結(jié)果與試驗(yàn)結(jié)果的偏差越大，但是Tstep2與Tstep1的計算結(jié)果相差無幾，可以說當(dāng)葉輪每個時間步轉(zhuǎn)過角度小于3°時，計算結(jié)果基本達(dá)到了時間步長無關(guān)性。

圖6 葉片法向力和切向力系數(shù)隨位置角的變化(λ=5.0)

4 結(jié) 論

通過對比以上計算結(jié)果可以得到以下結(jié)論:

(1)通過與試驗(yàn)結(jié)果對比，驗(yàn)證了CFD方法在垂直軸水輪機(jī)性能預(yù)報方面應(yīng)用的可靠性；

(2)湍流模型選擇方面，和其他湍流模型相比，SST湍流模型計算精度更高；

(3)當(dāng)滿足時基本滿足網(wǎng)格無關(guān)性要求；

(4)時間步長設(shè)定時，當(dāng)葉輪每個時間步轉(zhuǎn)過角度小于3°時，計算結(jié)果基本達(dá)到了時間步長無關(guān)性要求。

[1] Horiuchi, K., Ushiyama, I., and Seki, K., Straight Wing Vertical Axis Wind Turbine: A Flow Analysis, Wind Engineering, Vol. 29, No.3, pp.243-252, May 2005.

[2] Guerri, O., Sakout, A. and Bouhadef, K., Simulations of the Fluid Flow around a rotating Vertical Axis Wind Turbine, Wind Engineering Vol.31, no. 3, pp 149-163, 2007.

[3] 孫 科, 豎軸H型葉輪及導(dǎo)流罩流體動力性能數(shù)值模擬[J].哈爾濱工程大學(xué),2008.

[4] 姜子豪，李浩，闞江明.減振降噪在機(jī)械設(shè)計中的應(yīng)用[J].森林工程，2013，29(5)：55-59.

[5] Hwang, S., Lee, Y. H., and Kim, S. J., Optimization of Cycloidal Water Turbine and the Performance Improvement by Individual Blade Control, Applied Energy, Vol. 86, Issue 9, 2009.

[6] Ansys Inc: “ANSYS CFX User Manual”, 2005

[7] Deglaire, P., Engbom, S., Agren, O., and Bernhoff, H., Analytical Solutions for a Single Blade in Vertical Axis Turbine Motion in Two-Dimensions, European J. of Mechanics B/Fluids, Volume 28, Issue 4, pp. 506-520, July-August 2009.

[8] Tchon, Ko-Foa and Paraschivoiu, I., Navier-Stokes Simulation of the Flow Around an Airfoil in Darrieus Motion, Journal of Fluids Engineering, Vol. 116, pp. 870-876., ASME, 1994.

[9] 汪魯兵.豎軸潮流水輪機(jī)水動力性能理論與實(shí)驗(yàn)研究[J]. 哈爾濱工程大學(xué),2006.