韓端鋒，張 媛，高良田，孫 帥

(哈爾濱工程大學(xué)船舶工程學(xué)院，黑龍江哈爾濱150001)

0 引言

隨著CFD技術(shù)的迅速發(fā)展，利用CFD對(duì)水輪機(jī)進(jìn)行水動(dòng)力性能數(shù)值模擬、葉片優(yōu)化等方面的研究日趨成熟，計(jì)算精度也越來(lái)越高。

潮流能水輪機(jī)水動(dòng)力性能直接關(guān)系到水輪機(jī)發(fā)電功率的大小，提高水輪機(jī)的水動(dòng)力性能對(duì)水輪機(jī)的工程應(yīng)用具有重要意義。試驗(yàn)分析、CFD方法和理論分析3種方法是目前研究潮流能水輪機(jī)性能的主要方法。英國(guó)Lunar公司［1］運(yùn)用CFD方法分析了導(dǎo)流罩水輪機(jī)的水動(dòng)力性能，對(duì)導(dǎo)流罩水輪機(jī)進(jìn)行了1∶20水槽模型試驗(yàn)，驗(yàn)證了CFD方法的準(zhǔn)確性;英國(guó)南安普頓大學(xué)可再生能源實(shí)驗(yàn)室的Bahaj等人從獲能原理和試驗(yàn)分析的角度全面研究了水輪機(jī)特性，系統(tǒng)研究了水輪機(jī)效率、空化現(xiàn)象及波流耦合狀態(tài)下水輪機(jī)水動(dòng)力學(xué)性能，Bahaj［2］等的研究為潮流能水平軸水輪機(jī)的研究奠定了理論和試驗(yàn)基礎(chǔ)。國(guó)內(nèi)方面，北京風(fēng)能分析小組和Glasgow海洋能小組［3］運(yùn)用三維CFD模型研究了陣列水輪機(jī)性能與排布，分析了不同排布方式水輪機(jī)的水輪機(jī)水動(dòng)力性能，描繪了機(jī)組尾流場(chǎng)的速度場(chǎng)。本文通過(guò)CFD方法，研究水輪機(jī)對(duì)周?chē)鲌?chǎng)的影響，主要為一種新型的海洋新能源綜合發(fā)電平臺(tái) (見(jiàn)圖1)的總布置提供必要的參考。

圖1 海洋新能源綜合發(fā)電平臺(tái)Fig.1 Marine new energy integrated power generation platform

1 數(shù)學(xué)模型

計(jì)算流體動(dòng)力學(xué)是通過(guò)計(jì)算機(jī)數(shù)值模擬及顯示圖像對(duì)包含有流體流動(dòng)等物理現(xiàn)象的系統(tǒng)做分析。其主要思想是用一系列有限個(gè)離散點(diǎn)上的變量值的集合來(lái)代替原來(lái)在空間及時(shí)間域上連續(xù)的物理量的場(chǎng)，通過(guò)一定的方式和原則建立這些離散點(diǎn)上場(chǎng)變量之間關(guān)系的方程組，然后通過(guò)求解這些方程組得到場(chǎng)變量的近似值。它能夠不受實(shí)際試驗(yàn)所需各種條件的限制，對(duì)計(jì)算對(duì)象進(jìn)行離散和分析得出變量的近似值。本文主要是針對(duì)水輪機(jī)的流場(chǎng)進(jìn)行分析。

1.1 控制方程

流體動(dòng)力學(xué)的基本控制方程包括:質(zhì)量守恒方程［3］(連續(xù)性方程)、動(dòng)量守恒方程和能量守恒方程。本文中的計(jì)算首先假定流體為不可壓縮的流體，不涉及到求解能量守恒方程，對(duì)于牛頓流體而言，質(zhì)量守恒方程和動(dòng)量守恒方程的表達(dá)式為:

連續(xù)方程:

動(dòng)量方程:

其他變量輸運(yùn)方程:

這里的(i，j)指標(biāo)取值范圍是 (1，2，3)。根據(jù)張量的有關(guān)規(guī)定，當(dāng)某個(gè)表達(dá)式中一個(gè)指標(biāo)重復(fù)出現(xiàn)2次，則表示要把該指標(biāo)的取值范圍內(nèi)遍歷求和。

這里，τij實(shí)際對(duì)應(yīng)了6個(gè)不同的雷諾應(yīng)力項(xiàng)，即3個(gè)正應(yīng)力和3個(gè)切應(yīng)力。

由式(1)～式(3)構(gòu)成的方程共有5個(gè)方程 (雷諾應(yīng)力為3個(gè)方程)，而6個(gè)雷諾應(yīng)力再加上原來(lái)的5個(gè)時(shí)均未知量 (ui，p，φ)，總共有11個(gè)未知量，因此方程不封閉，必須引入新的湍流模型 (方程)，才能使方程組封閉。這就涉及到后面的湍流模型。

1.2 湍流模型

湍流是指雷諾數(shù)大于臨界值時(shí)，流動(dòng)特征發(fā)生變化，流動(dòng)成無(wú)序的這種狀態(tài)。它是一種高度非線性的流動(dòng)。為了研究這種流動(dòng)，通常把它看做時(shí)間平均流動(dòng)和瞬時(shí)脈動(dòng)流動(dòng)這2種流動(dòng)的疊加。引入雷諾平均法經(jīng)過(guò)轉(zhuǎn)換就得到式(3)的雷諾時(shí)均N－S方程 (RNAS)，然后通過(guò)引入湍流模型來(lái)封閉這個(gè)方程。

標(biāo)準(zhǔn)k－ω模型由于會(huì)在固壁上出現(xiàn)非物理的奇點(diǎn)，且不易處理避免附近變化劇烈的情況，Wilcox的k－ω被引入用于分析剪切流以及低雷諾數(shù)、可壓縮性的流體。其控制方程為:

式中:Gk與Gω分別為層流速度梯度與湍流頻率產(chǎn)生的湍動(dòng)能;Yk與Yω為擴(kuò)散產(chǎn)生的湍流;Sk與Sω為用戶定義的源項(xiàng)。

該模型對(duì)于低雷諾數(shù)下的流動(dòng)比標(biāo)準(zhǔn)k－ε模型具有更高的精度和數(shù)值穩(wěn)定性。

1.3 數(shù)值計(jì)算方法

采用有限體積法對(duì)三維不可壓縮RANS方程進(jìn)行離散，壓力速度耦合采用SIMPLEC算法進(jìn)行計(jì)算。壓力、密度以及動(dòng)量各項(xiàng)采用一階格式，為保證計(jì)算的收斂性，亞松弛因子均設(shè)置為0.2，設(shè)置殘差小于0.000 01時(shí)結(jié)果收斂。

2 計(jì)算模型及網(wǎng)格劃分

2.1 水輪機(jī)模型的建立和計(jì)算域的劃分

數(shù)值模擬計(jì)算的第一步是建立精確的幾何模型。它需要遵循在復(fù)合實(shí)際情況的前提下盡可能的簡(jiǎn)化模型。本文利用Catia軟件，根據(jù)表1中提供的數(shù)據(jù)，建模得到水輪機(jī)的幾何模型，如圖2(a)所示?？紤]到水輪機(jī)的旋轉(zhuǎn)運(yùn)動(dòng)，計(jì)算域分為2個(gè)區(qū)域，內(nèi)域和外域，內(nèi)域?yàn)樗啓C(jī)旋轉(zhuǎn)區(qū)域，內(nèi)含水輪機(jī)，外域?yàn)楣潭▍^(qū)域，為消除外域邊界對(duì)水輪機(jī)水動(dòng)力性能計(jì)算的影響，根據(jù)計(jì)算要求及模型尺寸建立恰當(dāng)?shù)目刂朴?，使其在滿足周?chē)鲌?chǎng)影響可忽略的情況下讓控制域盡可能的小。這樣可以在保證計(jì)算結(jié)果準(zhǔn)確的前提下盡可能的提高計(jì)算效率。本文取葉片和槳轂的表面作為內(nèi)邊界面，取直徑約為水輪機(jī)直徑4倍的圓柱體表面看作是無(wú)窮遠(yuǎn)邊界，作為外邊界面，取水輪機(jī)直徑1.1倍的小圓柱作為小域，便于進(jìn)行局部加密。具體模型如圖2(b)所示。

表1 潮流水輪機(jī)模型主要參數(shù)Tab.1 Geometry coefficients of the turbine

圖2 水輪機(jī)幾何模型和計(jì)算域Fig.2 Geometry model and calculation domain of turbine

2.2 網(wǎng)格劃分

網(wǎng)格的離散劃分是數(shù)值模擬中非常重要的工作，其好壞的程度將直接影響到求解的效率和精度。網(wǎng)格又分為結(jié)構(gòu)和非結(jié)構(gòu)網(wǎng)格兩大類。結(jié)構(gòu)網(wǎng)格適用于簡(jiǎn)單的模型，而對(duì)于復(fù)雜模型一般采用非結(jié)構(gòu)網(wǎng)格。所以對(duì)于流場(chǎng)計(jì)算問(wèn)題一般采用非結(jié)構(gòu)網(wǎng)格，但是，全部使用非結(jié)構(gòu)網(wǎng)格不僅網(wǎng)格質(zhì)量不是很好，大扭曲率的網(wǎng)格數(shù)較多，影響計(jì)算精度，也影響計(jì)算效率。對(duì)于這種情況，一般把整個(gè)流場(chǎng)劃分成幾個(gè)子域，然后根據(jù)具體情況在每個(gè)子域上選則適當(dāng)?shù)木W(wǎng)格類型。此外，網(wǎng)格的數(shù)目也是一個(gè)問(wèn)題，網(wǎng)格過(guò)梳會(huì)導(dǎo)致結(jié)果不精確甚至錯(cuò)誤;網(wǎng)格過(guò)密則會(huì)使計(jì)算時(shí)間過(guò)長(zhǎng)，計(jì)算難于收斂。因此在網(wǎng)格劃分過(guò)程中，我們需要合理安排網(wǎng)格數(shù)量，既要避免流場(chǎng)變化緩和區(qū)域的計(jì)算資源浪費(fèi)又要避免流場(chǎng)變化大的區(qū)域網(wǎng)格過(guò)疏造成計(jì)算不準(zhǔn)確，確保在節(jié)點(diǎn)數(shù)一定的情況下使用的網(wǎng)格模型能最大限度提高計(jì)算精度。另外，網(wǎng)格的劃分是一項(xiàng)非常繁瑣耗費(fèi)時(shí)間的工作，需要不斷的嘗試，直到生成復(fù)合要求的網(wǎng)格。本文采用混合式網(wǎng)格，在小域內(nèi)使用非結(jié)構(gòu)網(wǎng)格，并對(duì)水輪機(jī)葉面及槳轂進(jìn)行加密處理，在大域使用結(jié)構(gòu)網(wǎng)格，提高計(jì)算效率，縮短計(jì)算時(shí)間，具體網(wǎng)格如圖3所示。

圖3 網(wǎng)格的劃分Fig.3 Mesh of the turbine

2.3 邊界條件及參數(shù)的設(shè)定

邊界條件是指位于求解域邊界上的位置變量或其一階導(dǎo)數(shù)隨時(shí)間或者位置變化的規(guī)律只有選定合適的邊界條件，才能獲得流場(chǎng)正確的解。在Fluent中的基本邊界條件有流動(dòng)進(jìn)出口邊界，給定壓力邊界，壁面邊界，對(duì)稱邊界，周期性邊界。

進(jìn) (出)口邊界指在進(jìn) (出)口指定流動(dòng)參數(shù)的情況，一般可以設(shè)置速度、壓力、質(zhì)量進(jìn) (出)口邊界，針對(duì)可壓流動(dòng)的質(zhì)量進(jìn) (出)口邊界在這里不適用。由于流場(chǎng)進(jìn)行數(shù)值計(jì)算時(shí)，壓力有相對(duì)值表示，實(shí)際求解相對(duì)于進(jìn)口壓力而言所以有時(shí)須給定壓力參考值。壁面邊界作為最常用的邊界，采用它需修正壓力方程，尤其對(duì)于湍流計(jì)算，近壁面湍流演變?yōu)閷恿?，需?duì)離散方程的源項(xiàng)做特殊處理。對(duì)于周期性運(yùn)動(dòng)采用周期性邊界可以節(jié)約時(shí)間，提高計(jì)算效率。

由于水輪機(jī)在流體中運(yùn)動(dòng)時(shí)，流體不能與槳葉(轂)面分離又不能穿透槳葉 (轂)面，因此必須要求流體質(zhì)點(diǎn)在槳葉 (轂)面運(yùn)動(dòng)的法向分量與槳葉 (轂)面運(yùn)動(dòng)的法向分量相同。其次，水輪機(jī)置于無(wú)限大的流場(chǎng)中運(yùn)動(dòng)時(shí)，遠(yuǎn)處的流體基本不受擾動(dòng)，此時(shí)可認(rèn)為流體仍保持靜止?fàn)顟B(tài)。再次，由于不對(duì)尾流做研究可以假定尾渦面的厚度為零并且無(wú)流體和壓力降穿過(guò)尾渦面并且根據(jù)Kelwin渦量守恒條件得知流場(chǎng)中的渦不會(huì)消失。最后，要滿足庫(kù)塔條件即認(rèn)為隨邊葉面和葉背的壓力差為0，隨邊速度為有限值。

因此，針對(duì)此流場(chǎng)，在Fluent中將來(lái)流的入口設(shè)為速度入口，根據(jù)海域流速情況，選擇海域平均流速2.2 m/s為入口速度，將出口設(shè)為自由流動(dòng)。此外，將大域即整個(gè)流場(chǎng)設(shè)為無(wú)滑移壁面，而將水輪機(jī)設(shè)為滑移壁面與小域一起運(yùn)動(dòng)。

在設(shè)定邊界條件后，選擇SIMPLEC算法，剪切壓力傳輸 (SST)k－ω模型形成封閉的控制方程，然后通過(guò)設(shè)置求解的控制參數(shù)，初始化流場(chǎng)解得數(shù)值模擬的計(jì)算結(jié)果。

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

3.1 水輪機(jī)影響域計(jì)算

圖3為水輪機(jī)軸向截面流速分布云圖，該圖直觀地顯示了水輪機(jī)在尾流區(qū)域的流動(dòng)速度變化情況。由圖中可看出，水輪機(jī)對(duì)軸向方向上的速度影響范圍較大，在徑向方向上影響范圍較小;在葉梢位置形成一高速區(qū)，由于輪轂的阻礙作用，軸向速度迅速減小;尾流在水輪機(jī)的下游持續(xù)較長(zhǎng)，由圖5可知，距離輪轂中心位置越近，則速度恢復(fù)較慢;輪轂中心線處，由于輪轂的阻礙作用，軸向速度迅速減小為0，然后沿尾流方向，速度逐漸恢復(fù);在距離水輪機(jī)1.5D時(shí)，尾流速度為入口流速的63%;在距離水輪機(jī)2.5D時(shí)，尾流速度為入口流速的88%;葉尖位置處，在距離水輪機(jī)4D時(shí)，尾流速度為入口流速的93%;距離輪轂中心線0.5R處，沿尾流方向，速度先減小后恢復(fù)，其恢復(fù)較快。

圖4 單個(gè)水輪機(jī)影響域計(jì)算結(jié)果Fig.4 The influence domain calculation results of single turbine

圖5 不同位置軸向速度分布Fig.5 The distribution of axial velocity

3.2 水輪機(jī)間距對(duì)受力的影響

在相同的條件下，改變兩水輪機(jī)的橫向距離，進(jìn)行相關(guān)模擬計(jì)算［4］，得到水輪機(jī)距離分別為葉片直徑的1.5，2，3倍時(shí)水輪機(jī)受力情況，表2為兩水輪機(jī)受到的軸向力和轉(zhuǎn)矩值。

表2 不同間距下兩水輪機(jī)所受軸向力和轉(zhuǎn)矩值Tab.2 The axial force and the torque

由表2可看出，兩水輪機(jī)距離大于2倍直徑時(shí)，各水輪機(jī)的受力值在數(shù)值上接近單個(gè)水輪機(jī)的受力值，水輪機(jī)的距離對(duì)軸向力和轉(zhuǎn)矩的影響非常小，可忽略不計(jì);而當(dāng)水輪機(jī)間距為1.5 D時(shí)，對(duì)水輪機(jī)性能影響較大，兩水輪機(jī)軸向力及轉(zhuǎn)矩均有所減小。

3.3 水輪機(jī)間距對(duì)尾流場(chǎng)的影響

由圖6可看出，當(dāng)2個(gè)水輪機(jī)的間距為2倍直徑時(shí)，兩水輪機(jī)尾流速度之間幾乎沒(méi)有干擾;而當(dāng)水輪機(jī)間距縮小為1.5 D時(shí)，兩水輪機(jī)尾流方向速度存在較大的相互干擾，導(dǎo)致受影響的尾流區(qū)域相對(duì)變大。

為更好地描述水輪機(jī)工作的性能參數(shù)，選擇以下特征系數(shù)來(lái)表示:

水輪機(jī)從海流中獲取的能量用能量系數(shù)Cp表示，是評(píng)價(jià)水輪機(jī)性能最重要的參數(shù)，其表達(dá)式為:

式中:P為功率;V為水流速度;R為葉片半徑。

為研究水輪機(jī)在不同轉(zhuǎn)速下的特性，定義葉尖速比為λω，其表達(dá)式為:

式中ω為水輪機(jī)角速度。

圖6 水輪機(jī)間距對(duì)尾流場(chǎng)的影響Fig.6 The influence of turbine spacing on the wake flow field

圖7 水輪機(jī)性能曲線Fig.7 The performance of turbine

由圖7可知，采用CFD數(shù)值模擬得到的結(jié)果與實(shí)驗(yàn)值［6］趨勢(shì)相同，最大誤差在5%左右;當(dāng)尖速比為3.4時(shí)，試驗(yàn)值與數(shù)值計(jì)算值的Cp均開(kāi)始下降。CFD方法可以較為準(zhǔn)確得出水輪機(jī)的性能趨勢(shì)和最大效率，相比于試驗(yàn)方法，CFD方法更加簡(jiǎn)潔，成本較低，且能夠更加準(zhǔn)確地捕捉水輪機(jī)周?chē)鲌?chǎng)細(xì)節(jié)，對(duì)工程應(yīng)用具有一定的指導(dǎo)作用。當(dāng)2個(gè)水輪機(jī)間距為2D時(shí)，水輪機(jī)效率較單個(gè)水輪機(jī)效率略有減小，而當(dāng)間距減小為1.5D時(shí)，效率進(jìn)一步下降，尤其是當(dāng)尖速比增大時(shí)，效率下降明顯。說(shuō)明隨著水輪機(jī)轉(zhuǎn)速的增加，水輪機(jī)間相互干擾加強(qiáng)。

4 結(jié)語(yǔ)

本文首先對(duì)流體流動(dòng)的基本概念與相關(guān)流動(dòng)模型做了簡(jiǎn)單介紹，對(duì)Fluent軟件的基本結(jié)構(gòu)、數(shù)值求解方法和求解過(guò)程做了概述。其次，采用k－ε湍流模型結(jié)合混合網(wǎng)格技術(shù)開(kāi)展了水輪機(jī)敞水性能預(yù)報(bào)，計(jì)算了水平軸水輪機(jī)的水動(dòng)力性能，并對(duì)水輪機(jī)周?chē)鲌?chǎng)情況作了簡(jiǎn)要分析，討論了水輪機(jī)的安放間距對(duì)水輪機(jī)性能和流場(chǎng)的影響，研究發(fā)現(xiàn)，當(dāng)水輪機(jī)間距為1.5 D時(shí)，水輪機(jī)之間的相互影響較大，水輪機(jī)性能下降，間距為2 D時(shí)，水輪機(jī)相互之間幾乎沒(méi)有影響。

［1］ KERWIN JE，TAYLOR T E，BLACK SD，et al．A JSME Fluids Engineering，July 6－11，2003，Honolulu，Hawaii:ASME，2003:1－7．

［2］ BAHAJ A S，BATTEN W M J，MCCANN G．Experimental verifications of numerical predictions for the hydrodynamic performance of horizontal axis marine current turbines．Renewable Energy，2007，32:2479 －2490．

［3］ BAI L，ROY R G S，GREGORY D．Investagation of influence of array arrangement and spacing on tidal energy converter performance using a 3－Dimensional CFD model．

［4］ 陳存福．潮流能水平軸水輪機(jī)葉片優(yōu)化及水動(dòng)力性能研究［D］．青島:中國(guó)海洋大學(xué)，2012．CHEN Cun-fu．Study on blade optimization and hydrodynamic performance of horizontal axis marine current turbine［D］．Qingdao:Ocean University of China，2012．

［5］ 劉志華，熊鷹，葉金銘，等．基于多塊網(wǎng)格混合的RANS

方法預(yù)報(bào)螺旋槳水動(dòng)力性能的研究［J］．水動(dòng)力學(xué)研究與進(jìn)展，2007，22(4):450 －456．LIU Zhi-hua，XIONG Ying，YE Jin-ming，et al．Based on hybrid RANSmethod of forecast of the grid block propeller hydrodynamic performance research［J］．Water Dynamics Research and Development，2007，22(4):450 －450．

［6］ 安佰娜．潮流能發(fā)電場(chǎng)尾流場(chǎng)數(shù)值模擬及其多機(jī)組影響規(guī)律研究［D］．青島:中國(guó)海洋大學(xué)，2012．An Bai-na．Numerical Analysis of the Wake and Interaction Effect Analysis for tidal arrays in the Tidal Farm［D］．Qingdao:Ocean University of China，2012．

［7］ 王福軍．計(jì)算流體動(dòng)力學(xué)分析［M］．北京:清華大學(xué)出版社，2004．WANG Fu-jun．Computational fluid dynamics analysis［M］．Beijing:Tsinghua University Press，2004．