田琳琳，趙 寧，鐘 偉，胡 偶

(南京航空航天大學航空宇航學院，江蘇南京 210016)

0 引言

風力機是將自然界的風能轉化為機械能并獲得電能的裝置。由于能量的轉移，風經過旋轉的風力機之后，流動情況發(fā)生了很大的變化:風速減小、湍流強度增加、出現(xiàn)了明顯的風剪切層。風速減小會使下游風力機的輸出功率降低，尾跡附加的風剪切和強湍流會影響下游風力機的疲勞載荷、使用壽命和結構性能。經過一段距離之后，在周圍氣流的作用下，風速逐漸得到恢復。這就是風力機的尾流效應［1］。為了使有限面積的風場產生最大的經濟效益，風電場開發(fā)商盡可能多在風場內布置風力機，這樣上游機組會對下游機組產生較大的尾流效應，導致下游機組發(fā)電量減少，進而降低整個風場的總輸出功率。因此，開展尾流效應對風力機性能影響的研究對合理布置風力機，減小尾流對下游的影響，提高風能資源的利用效率，減小風電場的占地面積，最終使風電場的經濟性達到最佳起著重要的作用。

風力機的尾流區(qū)分為近尾流區(qū)和遠尾流區(qū)。近尾流區(qū)的研究著眼于功率提取的物理過程和風機性能，而遠尾流區(qū)的研究重點是尾流模型、地形影響、湍流模型的恰當選取，更著重于研究風電場機組間的相互干擾［1］。

許多專家學者對風力機的尾流效應開展了廣泛的計算研究?？梢詫⑦@些研究大致分為兩類:一類是尾流模型研究，這些模型是由學者提出的簡化的尾流模型，然后利用實驗數據檢驗模型并進行修正，主要用來解決工程問題。如風資源分析軟件WAsP采用的Park模型［2］以及Ainslie提出的渦粘性尾流模型［3］等等。另一類是流動機理研究，該部分主要是回歸基礎研究，旨在揭示流動的內在規(guī)律和本質。如早期的動量－葉素理論［4］和渦尾跡方法［5］。之后隨著計算機的發(fā)展，又出現(xiàn)了很多基于CFD的方法，例如Ivanell等學者［6］使用了制動線(ACL)方法研究了單個風力機的近尾跡區(qū)的流動狀況，Sankar等學者構造了N－S方程、全勢能方程與自由渦尾跡相混合的方法用于求解水平軸風力機周圍的非定常粘性流場［7］。S?rense［8］等使用 RANS 方程進行流場的全三維數值模擬。

本文提出了兩種方法用來研究單個風力機的遠尾跡特征。第一種方法為修正的Park模型，由初始的一維Park模型進行三角函數修正得到。Park模型由于計算量小且相對準確［9］，在工程上得到了廣泛的應用，經過修正的Park模型計算量也很小，且在徑向方向的速度型更接近于真實的流場分布情況。第二種方法是結合制動盤理論和CFD［10］，將風力機簡化成為一個無限薄的產生壓力躍升的制動盤，采用FLUENT商用軟件提供的Fan邊界對流場進行數值模擬。該方法較之建立真實的風力機模型并進行模擬來說，降低了計算量和計算復雜度，為風力機的微觀選址提供很好的理論支持。最后將本文的兩種數值方法的結果與風洞實驗結果以及Park模型的預估結果進行對比分析，發(fā)現(xiàn)這兩種方法都能夠很好地預測風力機的尾流流場。

1 研究模型

1989年 Garrad Hassan在 Marchwood工程實驗室(MEL)開展了一系列的風力機尾流的實驗研究［11］，所選用的實驗對象為模型縮比尺度為1/160的水平軸風力機。對于真實尺度的風力機來說，能反映其性能的重要的參數就是推力系數，所以在風洞試驗中，觀察了風力機在三個不同尖速比下的氣動性能。三個不同尖速比 λ =2.9，λ =4.0 ，λ =5.1 對應的軸向推力系數分別為 CT=0.62，CT=0.79 ，CT=0.85 。

本文選用的計算模型為Hassan實驗的縮比模型所對應的真實模型，即風輪直徑為43.2m，輪轂高度為50m，地表粗糙度為0.075m，相當于風力機樹立于短的農作物覆蓋的真實地形之上。

2 計算方法

2.1 Park模型的修正

Park模型是由學者Jensen于1986年提出的一維尾流模型［2］，該模型假設在風力機下游，尾流是線性膨脹的，Park模型的原理圖如圖1所示。根據質量守恒定理推導出在風力機下游x位置處速度為:

其中u0為來流風速，a為軸流誘導因子，與推力系數CT有關，r1為下游x位置處尾流區(qū)域的半徑，α為尾流擴散系數，表示尾流的膨脹速率。a、r1和α的公式如下:

其中，rr為風輪半徑，z為風力機輪轂高度，z0為地表粗糙度。

圖1 尾流模型原理圖Fig.1 Schematic of wake model

大量的實驗結果以及數值模擬表明，在真實的尾流流場中，沿風力機徑向的速度分布應該為拋物型的［12］，這與三角函數的圖形有著很大的相似性，故在此使用三角函數對一維的Park模型進行修正，使其為二維模型，保證修正后的徑向速度型為拋物型。修正后的Park模型示意圖如圖1所示。

令y方向為徑向方向，假設修正后的速度分布符合標準的三角函數關系式:

下面將分別求解A，k以及b。根據余弦函數周期的定義可得:

假設余弦函數取最大值時所對應的徑向位置為r1即膨脹尾流的半徑，則有:

基于修正前后動量守恒的原則，也就是說修正后的總動量與Park模型推導出的總動量相等，即速度型下所包含的面積相等，有:

聯(lián)立(6)～(8)式分別得到:

可得u*的表達式為:

所以經過修正的Park模型分為兩步:

2.2 制動盤理論結合CFD方法

當空氣經過風輪盤時，由于風力機的阻塞作用，使得流管內的風速降低，根據質量守恒原理，流管的橫截面積就會膨脹以適應減速的空氣，我們將風輪所在平面稱為制動盤［12］。

由經典的Rankine－Froude原理可知，風經過制動盤時，會在制動盤前后產生壓差，壓差的表達式如下:

其中為△P空氣經過制動盤時產生的壓力降，ρ為空氣密度，U0為來流速度5.3m/s，a為軸流誘導因子，見式(2)。

把風力機簡化為產生壓力降的制動盤，這與商用軟件FLUENT提供的Fan的物理模型具有相似性，所以可借助Fan邊界來模擬風輪，計算Fan對整個流場的影響。在FLUNT中，把Fan假設為一個無限薄的面，流體經過這個面時出現(xiàn)壓力躍降，而躍降的大小是速度的多項式函數，壓強躍降函數的表達式為:

在式(13)中△P為壓強躍降量，fn為多項式系數，v為速度。這樣做的目的是把Fan前后的速度變化轉化成壓差項反映到流動方程中。

在本文的計算中，將風力機簡稱為一個無限薄的制動盤，并根據風力機的推力系數與風速關系的特征曲線設定制動盤前后的壓差。然后在全流場內求解帶湍流模型的N－S方程。

3 計算結果與分析

3.1 風輪中心軸向速度分布

能否準確地計算尾流軸向動量虧損是判斷一個模型是否能夠很好地預測尾流發(fā)展情況的一個重要指標。圖2給出了在三個不同尖速比下，制動盤結合CFD方法、Park模型以及風洞實驗這三種方法分別得到的沿風輪中心的速度輪廓圖。在圖中，x方向是來流速度方向，D為風輪直徑，U是尾流的當地速度，U0是無窮遠處的來流速度。由于Park模型在風輪下游的兩倍風輪直徑之后才有效，所以本文選擇從風力機下游的2D距離之后對三種方法進行分析比較。圖中很明顯地反映了風力機的下游速度會下降，經過一段距離之后，逐步恢復到來流速度。且隨著尖速比的增加，推力系數就會增加，尾流效應就越明顯。在下游6D之前，Park模型和本文的制動盤結合CFD方法都偏大很多，是因為這兩種方法都沒有考慮到風力機的旋轉效應，而真實的風洞試驗中，風力機是旋轉的，旋轉會引起風力機的下游產生很強的湍流以及尾渦，尾渦通過發(fā)展、破碎等過程會耗散掉空氣中很多的能量，導致下游的速度進一步降低。Helmis等學者［13］指出風洞實驗過高預估了近尾跡區(qū)的尾流效應，使分析得出的風輪中心的速度虧損較大。這也是本文方法的結果與風洞試驗結果在近尾流區(qū)相差較大的原因之一。

開展尾流效應研究的目的是在風電場內合理布置風力機，風力機的安裝間距要滿足風場總體效益最大化的目標。通過對國內外風場多年的建設經驗分析，風力發(fā)電機組安裝間距在盛行風向上選擇為6～10倍風輪直徑，在垂直盛行風向上選擇為3～5倍的風輪直徑較為合適。在下游距離為6～15倍風輪直徑之間，本文的制動盤結合CFD方法的計算結果與實驗值吻合的很好，而Park模型的計算結果偏高。這說明本文的計算方法可以為風電場的微觀選址提供很好的理論依據。

圖2 風輪中心軸向速度分布Fig.2 Velocity deficit on the wake centerline

在風力機下游15倍風輪直徑之后，Park模型的計算值比實驗值偏大，本文的計算結果偏小，說明尾流還沒有來得及完全恢復。

3.2 下游位置的徑向速度分布

圖3 尖速比2.9時，風力機下游的不同位置處的徑向速度分布Fig.3 The crosswind profiles of mean velocity at hub height for the tip speed ratio of 2.9

圖3和圖4分別給出了不同尖速比下風力機下游的不同橫截面在風輪高度處的徑向速度分布。該圖明顯的顯示了尾流效應，且尖速比越大，尾流效應越明顯。從整體來看，Park模型以及本文的制動盤結合CFD方法的計算結果都大于實驗結果，這是因為在風力機下游的一段距離之后，壓差不再在整個流動中起主導作用，而是湍流強度以及尾渦的變化主導著整個流場的能量變化，在真實情況中，風力機下游的尾渦會逐漸耗散掉，導致能量的降低，進而速度會降低。在越靠近下游的地方，制動盤結合CFD方法以及Park模型結果會更接近風洞試驗值，尤其是在尖速比為2.9時，風力機下游10倍風輪直徑的地方，本文結果與風洞試驗值幾乎吻合。從圖中還可以看出經過修正的Park模型與風洞試驗結果和制動盤理論結合CDF方法計算結果吻合的較好，尤其是在風輪中心兩側，與實驗值更加接近。由于其計算量小，計算復雜度低，因此在今后可作為一項計算工具應用于工程實際。圖3和圖4還反映了制動盤結合CFD方法計算處的尾流寬度略大于Park模型的尾流寬度，更接近于實驗值。一個比較有趣的現(xiàn)象是在流管邊緣的地方，制動盤結合CFD方法計算的結果比實驗值和Park模型的值偏大2%，風出現(xiàn)了稍微的加速現(xiàn)象，這是因為在制動盤的周圍拖出來的尾渦從周圍流體中吸收了能量，導致速度的稍微上升。

圖4 尖速比5.1時，風力機下游的不同位置處的徑向速度分布Fig 4 The crosswind profiles of mean velocity at hub height for the tip speed ratio of 5.1

4 總結

本文根據GH公司提供的風力機數據，通過兩種方法對單個風力機的尾流進行了計算模擬。一種方法是構造新的尾流模型計算流場速度分布，另一種方法是結合制動盤理論和CFD，在全流場進行數值模擬。當風吹過風力機時，風速會下降，經過下游一段距離慢慢恢復到來流速度，且尖速比越大，尾流效應越明顯。從風輪中心軸向速度分布圖來看，在下游6倍風輪直徑之后，本文結果與風洞試驗值的最大誤差約為4%。徑向速度分布顯示了尾流的膨脹系數大概為0.08，越靠近下游，本文計算方法的精度就越高。雖然Park模型沒有精確預測尾渦的形狀，但是它近似給出了尾渦的發(fā)展過程，而本文構造的新尾流計算模型與風洞實驗值吻合較好。

將本文兩種計算方法的結果與風洞實驗以及Park尾流模型進行分析比較，顯示了本文能夠定性反映風力機周圍的物理流動現(xiàn)象。以上結果表明，本文的兩種計算方法適用于風力機的流場模擬，可作為工程應用的工具，用于風電場的流場模擬，為風力機的微觀選址提供一定的參考依據。

［1］VERMEER L J，S?RENSEN J N，CRESPO A.Wind turbine wake aerodynamics［J］.Progress in Aerospace(Sciences)，2003，39(6 －7):467 －510.

［2］JENSEN N O.A note on wind generator interaction［R］//Ris?National Laboratory，1983.

［3］AINSLIE J F.Development of an eddy viscosity model for wind turbine wakes［A］.Proceeding of 7 th BWEA Wind Energy Conference［C］，Oxford，1985.

［4］GLAUERT H.Airplane propellers［A］.In Aerodynamics theory.Dover Publication Inc:New York［M］，1963，4(Division L):169－360.

［5］BAREISS R，GUIDATI G，WAGNER S.Wake simulation for wind turbines with a free，prescribed and hybrid wake method［A］.Proc.10th IEA Symp.On the Aerodynamics of Wind Turbines［C］，Edinburgh，1996:23－28.

［6］IVANELL S，MIKKELSEN R，S?RENSEN J，HENNINGSON D.A stability analysis of the tip vortices of a wind turbine［J］.Wind Energy，2009.

［7］XU G，SANKER L N.Computational study of horizontal axis wind turbines［R］.AIAA Paper 99－0042，1999.

［8］GLAUERT H.Airplane propellers［A］.In Aerodynamics theory［M］.Dover Publication Inc:New York，1963，4(Division L):169－360.

［9］BARTHELMIE R，F(xiàn)OLKERTS，LARSEN G C，et al.Comparison of wake model simulation with offshore wind turbine wake profiles measured by sodar［J］.Journal of Atmospheric and Oceanic Technology，23(7):888-901.

［10］DZUIBINSKI Z，WIENZYSLAW S.Vortex ring state simulation using actuator disc［EB/OL］.http://www.scs－europe.net/conf/ecms2007/ecms2007－cd/ecms2007/ecms2007%20pdf/ese_0121.pdf.

［11］SCHLEZ W，TINDAL A，QUARTON D.GH windfarmer validation report［R］.Bristol:Garrad Hassan and Partners Ltd，2003.

［12］EBERLE W R.Wind flow characteristics in the wake of large wind turbines，Vol I，analytical development［R］.DOE/NASA CR－165380.

［13］賀德馨，等.風工程與工業(yè)空氣動力學［M］.北京:國防工業(yè)出版社，2006.

［14］HELMIS C G，PAPADOPOULOS K H，ASIMAKOPOULOS D N，PAPAGEORGAS P G，SOILEMES A T.An experimental study of the near wake structure of a wind turbine operating over complex terrain［J］.Solar Energy，1995，54(6):413-28.