賀廣零,佟 輝
(1.中國電力工程顧問集團華北電力設計院工程有限公司,北京市100120;2.大慶油田電力集團宏偉熱電廠,黑龍江省大慶市163411)
基于常規的紊流風速功率譜(如Von Karman 譜,Kaimal 譜,Davenport 譜等)很難準確預測風力發電機組的疲勞荷載,實際的疲勞荷載往往遠比預測值大。此外,也無法成功建立實測風速和風力發電機組動力響應(尤其是塔筒基底彎矩)之間的聯系。現場實測表明,旋轉葉片上的紊流風譜與固定點的紊流風譜大相徑庭,其能量分布發生了根本變化。總體上,是由于葉片旋轉效應導致作用在葉片上的風場結構發生變化所致。因此,為準確描述作用在旋轉葉片上的風荷載,必須考慮葉片旋轉效應。1955年,Rosenbrock[1]最早提出了考慮葉片旋轉效應的紊流風速模型。半個世紀以來,考慮葉片旋轉效應的紊流風速模型主要有:PNL 模型和SNL 模型。
為了能夠描述作用在旋轉葉片上的風速,Connell[2-3]將歐拉坐標系轉換到旋轉坐標系,并在該坐標系下建立了一種風速功率譜數學模型,即PNL模型。之后,Powell & Connell[4-7]、Connell[8]對PNL模型進行了持續完善,并推廣應用。然而,PNL 模型自身存在如下缺陷:(1)僅能模擬葉片上一點的風速,無法考慮風輪平面上不同點之間的相關性;(2)未能考慮葉片數量的影響;(3)需要源譜在進行Fourier 變換后能獲得解析表達式,要求過于苛刻。
與此同時,Veers[9-10]對風輪平面上各取樣點進行風場模擬以獲取風速時程,然后隨著葉片旋轉移動逐點提取瞬時風速,所提取的瞬時風速依序構成一組新的風速時程。與PNL 模型相比,SNL 模型簡單、直觀,能對整個風場進行有效模擬,且對源譜和相干函數沒有特殊要求,故而更具開放性。然而,該模型同樣存在一些弊端:(1)只是在風速層次上考慮了葉片旋轉效應,并未真正構建旋轉樣本譜;(2)進行風場模擬時,相位信息不太容易精確控制,容易導致較大的誤差;(3)葉片上點暫時不經過的位置也會進行風場模擬,計算效率較低。PNL 模型和SNL 模型具有共同的本質,是從不同角度來詮釋葉片的旋轉效應。簡單來說,PNL 模型是從Lagrange 描述的角度研究葉片旋轉效應,而SNL 模型則立足于Euler 描述。有關葉片旋轉效應,盡管已有學者進行了較長時間的探索,但其物理機制并非如人們所期望的那樣清晰,其物理意義及工程價值還有待進一步探索。
本文從一個新的角度來探索葉片旋轉效應的物理機制,推導基于物理機制的葉片旋轉Fourier 譜。在推導過程中,擯棄平穩過程概念和各態歷經假定,反映真實動力過程的特征。進一步地,在完成旋轉Fourier 譜參數分析的基礎上,明確其物理意義。
近年來,相當一部分學者認為[11-12]:風場中大部分漩渦尺寸要小于葉片尺寸,當該漩渦經過葉輪平面時,葉片將切割漩渦若干次,即所謂的"漩渦切割"現象(圖1)。換句話說,相對固定點的風速而言,旋轉葉片上的風速經歷了紊流風場的空間結構變化。因為葉片切割漩渦若干次,葉片將經歷相應次數紊流風場結構的變化,從而導致紊流風譜出現能量重分布。顯然,這種對葉片旋轉效應物理機制的闡述晦澀難懂,且難以用數學公式來衡量。

圖1 漩渦切割現象Fig.1 Eddy slicing phenomenon
為撥開葉片旋轉效應的重重迷霧,筆者將從另一角度探索其物理機制。一般來說,紊流風場是一個時變隨機場,即空間中任意一點的紊流風速不僅與該點的空間位置有關,而且還與時間有關。此外,風力發電機組紊流風場還具有其特殊性:葉片上任意一點的空間位置隨著葉片旋轉而不斷變化,從而使得紊流風場空間變化性更為復雜。為了能較好地刻畫作用在旋轉葉片上的風速時程,首先需要轉換思考角度,將坐標系切換至旋轉坐標系。為了能跟蹤旋轉葉片的軌跡,可有規律地沿著旋轉葉片的軌跡選取樣本點(圖2)。在旋轉葉片經過某樣本點時,提取該樣本點風速時程在該時刻的瞬時風速,依此類推,將提取的瞬時風速按既有的時間順序排列可構建一組新的風速時程。顯然,這組風速時程的共同特征在于:自始至終作用于旋轉葉片上某一點,而這組風速時程對應的紊流風速譜模型即為旋轉Fourier 譜。

圖2 風輪平面上取樣點位置Fig.2 Sampling location on rotor plane
在旋轉Fourier 譜物理過程分析中,需要獲得樣本互相關函數與Fourier 互譜之間的關系。在此,首先對二者之間的關系進行推導。值得說明的是,以下推導均基于有限時長的樣本時程,擯棄了平穩過程概念和各態歷經假定,反映了真實動力過程的特征。
對于2個不同的樣本時程x1(t)、x2(t + τ),其Fourier 變換可分別定義為

其相應的樣本Fourier 密度譜可分別定義為

同時,時程x1(t)、x2( )t 的樣本互相關函數可定義為

利用式(1)、(2),可得

故

即

同理可得

顯然,Fourier 互譜Fx12(n)可以表示為

在另一方面,引入相干函數概念,Fourier 互譜也可以表示為

顯然,相干函數(γ d(τ))為

一般情況下,(γ d(τ)) 取決于兩點之間的距離d(τ)。
依據式(12),式(13)可以改寫為

由上式可知,樣本互相關函數與Fourier 互譜構成Fourier 變換對。
樣本時程x(t)自相關函數R11(τ) 可定義為

依據(8)式,當x1(t) = x2(t)時可得


即

由上式可知,樣本自相關函數與Fourier 自譜構成Fourier 變換對。
由旋轉Fourier 譜的物理機制可知,旋轉Fourier互譜與經典的紊流風譜互譜存在本質的不同,已然不再是自譜與相干函數的簡單乘積。已知葉片以轉速n0勻速旋轉,半徑ri上的一點在t 時刻的風速時程樣本為xi,半徑rj上的一點在t + τ 時刻的風速時程樣本為xj(圖3)。d(τ)為兩點間距離(圖3),其表達式為

式中φ 為相位因子。當i、j 兩點處于同一片葉片上時,φ = 0 ;當i、j 兩點處于不同葉片上時,φ 則為2π/Nb(Nb為葉片數目)的整數倍。對于三葉片風力發電機而言,φ 始終維持為2π/3 。
由(14)式可知,Fourier 互譜與互相關函數構成Fourier 變換對。對Fourier 互譜[13]進行逆Fourier 變換,可得兩點間的互相關函數為


圖3 旋轉葉片相對位置關系Fig.3 Relative position of rotational blade
根據旋轉Fourier 譜的物理機制,旋轉葉片上任意兩點在不同時刻的相關函數可用兩點間的互相關函數來代替,即

值得注意的是,Rij(τ)為i 點與轉動時間τ 后j點的互相關函數。對互旋轉相關函數進行Fourier 變換可得到旋轉Fourier 互譜ij(n)為

將式(20)、(21)代入式(22)中

與前述分析類似,將γ(d(τ),n')進行Fourier 級數展開

式中:n0為葉片轉動頻率;k'm(n)為相干函數γ(d(τ),n)的Fourier 展開系數。

將(24)式代入(25)式中,有

亦即


進而可得

式中,k'm(n - mn0)為旋轉模態,本質上為相干函數的Fourier 展開系數。eimφ為葉片標識項,當φ = 0 時,為同一片葉片上的互譜;當φ ≠0 ,為不同葉片上的互譜。令ri= rj且相位因子φ = 0 ,可得旋轉Fourier自譜ii(n)為

以典型的1.25 MW 三葉片變槳距風力發電機組為例進行分析。依據風機廠家提供的資料,該機型主要參數如下:葉輪直徑為64 m,轉速為30 r/min(0.5 Hz),輪轂高度為68 m,10 m 高處平均風速為15 m/s,地面粗糙度為0.02 m,湍流積分長度73.5 m。
圖4 給出了旋轉Fourier 譜與樣本Fourier 譜之間的比較。由圖4 可知,旋轉Fourier 譜出現能量重分布:由低頻向高頻大幅度轉移,并在旋轉頻率n0的整數倍處出現峰值。由式(29)可知,旋轉Fourier 譜ij(n)是由無限個樣本Fourier 譜Fii(n)疊加而成,而這些樣本Fourier 譜在疊加之前先進行了葉片旋轉頻率n0整數倍平移。如此,旋轉Fourier 譜能量重分布的本質一目了然。

圖4 旋轉Fourier 譜與樣本Fourier 譜比較Fig.4 Comparison between rotational &sampling Fourier spectrum
不難判斷,旋轉Fourier 譜與計算點半徑、葉片旋轉頻率密切相關。同時,旋轉Fourier 譜以隨機Fourier 譜為源譜[14],自然相應地以10 m 高平均風速、地面粗糙度為主要影響因素。以下將針對各影響因素逐個進行參數分析。
(1)計算半徑的影響。圖5 分別給出了計算半徑為30,20,10,0 m 時的旋轉Fourier 譜曲線。由圖5不難發現:計算點半徑r 對旋轉Fourier 譜的影響:當r =30 m 時,旋轉效應最為顯著,能量重分布最為突出;r =20 m、r =10 m 旋轉效應依次減弱;r =0 m 時,葉片旋轉效應徹底消失。總體上,隨著計算半徑的遞減,葉片旋轉效應逐漸減弱。

圖5 計算點半徑對旋轉Fourier 譜的影響Fig.5 Influence of computing radius on rotational Fourier spectrum
(2)旋轉頻率的影響。圖6 展現了葉片旋轉頻率n0對旋轉Fourier 譜的影響。對于大型風力發電機組而言,葉片旋轉頻率范圍大致為0.2 ~0.6 Hz,不妨取0.2,0.4,0.6 Hz 為代表進行分析。旋轉Fourier 譜對旋轉頻率高度敏感,旋轉頻率的些許變化就會引起旋轉Fourier 譜曲線較大的調整,主要體現在以下2個方面:1)隨著旋轉頻率的增加,旋轉Fourier 譜峰值快速增大,谷值快速減少,葉片旋轉效應趨于顯著;2)旋轉Fourier 譜的峰值出現在葉片旋轉頻率的整數倍處。因此,隨著旋轉頻率的增加,旋轉Fourier 譜峰值整體向右平移。值得指出的是,三者在低頻階段基本重合,只是在高頻階段開始分化。

圖6 旋轉頻率對旋轉Fourier 譜的影響Fig.6 Influence of rotational frequency on rotational Fourier spectrum
(3)平均風速的影響。圖7 展示了輪轂高度處平均風速對vh旋轉Fourier 譜的影響。一般地,風力發電工程以10 m 高度處平均風速為代表,而風電項目則習慣以輪轂高度處平均風速為代表,二者存在一定的換算關系。該風力發電機組的切入風速、切出風速分別為4、25 m/s,因此在二者中間選擇了3個代表性風速(8、16、24 m/s)進行分析。為了獲得比較明顯的計算結果,選定計算點距輪轂中心30 m。不難發現,當輪轂高度處平均風速較小時,旋轉Fourier譜峰值較大,譜谷值較小,譜曲線形狀細長、尖銳。隨著平均風速的逐漸增大,旋轉Fourier 譜峰值減少,譜谷值增大,譜曲線形狀趨于平緩、豐滿。此外,在低頻階段,三者峰值對應的頻率存在一定差異。

圖7 平均風速對旋轉Fourier 譜的影響Fig.7 Influence of mean wind speed on rotational Fourier spectrum
(4)地面粗糙度的影響。圖8 體現了地面粗糙度z0對旋轉Fourier 譜的影響。地面粗糙度z0與地形地貌緊密相關。在分析過程中,賦予z0= 0. 005,0.07,0.3,1.0,分別對應于海岸區、開闊地形、郊區地形、城市中心四種地形。由圖8 可知,隨著地面粗糙度z0增大,旋轉Fourier 譜值整體向上平移。其原因在于增大地面粗糙度z0會使得剪切流動速度變大,進而導致脈動風標準差增大,最終引起旋轉Fourier 譜能量增加,在圖形上則表現為譜值向上平移。
(5)自譜與互譜的比較。圖9 為旋轉Fourier 互譜與自譜的比較。以不同半徑r1=30 m 和r2=20 m上的2 點之間的互譜為例,同時與r1= r2=30 m,r1=r2=20 m的自譜進行比較。互譜的形狀與自譜的形狀大體相似。在低頻段,互譜值介于2個自譜值之間;而在高頻段,互譜值衰減非常迅速,以至于遠小于2個自譜值。

圖8 地面粗糙度對旋轉Fourier 譜的影響Fig.8 Influence of surface roughness on rotational Fourier spectrum

圖9 旋轉Fourier 互譜與自譜的比較Fig.9 Comparison between cross &auto rotational Fourier spectrum
(1)旋轉Fourier 譜以隨機Fourier 譜為源譜,而隨機Fourier 譜本質上為物理隨機函數模型。物理隨機函數模型可以有效地避免經典現象學隨機過程模型的種種局限性[15],不僅可以給出風電機組隨機風場的完整數學描述,而且由于其具備明確的物理意義,使得對于風電機組隨機風場的實驗驗證成為可能,同時也有助于發現隨機性的客觀來源。
(2)PNL 模型和SNL 模型都是從時域角度來分析葉片旋轉效應,本文則從頻域的角度出發,建立了旋轉Fourier 譜物理模型。這不僅是思考角度的變化,更重要的是建立了旋轉Fourier 譜與其源譜之間的映射關系,這意味著旋轉Fourier 譜的源譜可以是任意一種合理的紊流風譜。同時,由于以隨機Fourier 譜為源譜,且該譜本質上為隨機函數,具有明確的物理機制,是迄今為止最為合理的紊流風速譜。當然,旋轉譜的源譜也可以是經典功率譜,如Karman譜、Kaimal 譜、Simiu 譜等,功率譜對應的旋轉譜可定義為旋轉樣本譜[16]。同時,若將來出現的更為合理的風速譜,則可直接納入本文提出的映射機制獲得更為精確的旋轉譜。因此,本文提出的旋轉譜不僅簡單、明確,而且具有很強的兼容性和普適性。
(3)旋轉Fourier 譜反映了風速作用在旋轉葉片上的物理機制,其本質上為作用在旋轉葉片上紊流風速時程對應的隨機Fourier 譜。旋轉Fourier 譜反映了作用在旋轉葉片上風速的時間、空間雙重變化,不僅體現了風速自身的脈動特性(時變特性),而且刻畫了葉片高度周期性變化引起的脈動風速調整(空間變化特性)。因此,旋轉Fourier 譜能更準確地預測風力發電機組極值荷載,尤其是疲勞荷載。
(4)旋轉Fourier 譜的突出優勢是:在保證分析結果可靠的前提下,能將葉片旋轉這一個運動學問題簡化為靜力學問題[13,17-18]。如何考慮葉片旋轉效應是風力發電機組風致隨機動力響應分析的核心問題。事實上,基于有限元法建模無法模擬葉片旋轉(剛體運動)。運動具有相對性,以旋轉葉片為參照物,視葉片風場在不停地反方向旋轉,基于這種思想構建旋轉Fourier 譜,顯然已經有效地考慮了葉片旋轉效應。因此,基于該譜進行風力發電機組風致動力響應分析時無須再次考慮該效應,這對風力發電機組結構分析而言,無疑是個很大的簡化。
基于風速作用于旋轉葉片的物理機制,對相干函數進行Fourier 級數展開,利用δ 函數的性質,推導了旋轉Fourier 自譜和互譜表達式。研究表明,旋轉Fourier 譜會出現能量重分布現象,且與計算半徑、葉片旋轉頻率、平均風速、地面粗糙度等因素密切相關。采用旋轉Fourier 譜,能夠更準確地預測風力發電機組極值荷載,尤其是疲勞荷載;在保證分析結果可靠的前提下,能將葉片旋轉這一個運動學問題簡化為靜力學問題。
[1]Rosenbrock H H. Vibration and stability problem in large turbines having hinged blades[D]. E.R.A.,Surrey,U.K.,1955.
[2]Connell J R. Turbulence Spectrum observed by a fast-rotating wind turbine blade [R]. Technical Report PNL-3426,Batelle Pacific Northwest Laboratory,Richland,WA,1980.
[3]Connell J R. The spectrum of wind speed fluctuations encountered by a rotating blade of a wind energy conversion system [R].Technical Report PNL-4083,Batelle Pacific Northwest Laboratory,Richland,WA,1982.
[4]Powell D C,Connell J R,George R L. Verification of theoretically computed spectra for a point rotating in a vertical plane[R].Technical Report PNL-5440,Battelle Pacific Northwest Laboratory,Richland,WA,1985.
[5]Powell D C,Connell J R. A model for simulation rotational data for wind turbine applications[R]. Technical Report PNL-5857,Battelle Pacific Northwest Laboratory,Richland,WA,1986.
[6]Powell D C,Connell J R. Review of wind simulation methods for horizontal-axis wind turbine analysis[R]. Technical Report PNL-5903,Battelle Pacific Northwest Laboratory,Richland,WA,1986.
[7]Powell D C,Connell J R. Verification of theoretically computed spectra for a point rotating in a vertical plane[J]. Solar Energy,1987,39(1):53-63.
[8]Connell J R. A primer of turbulence at the wind turbine rotor[J].Solar Energy,1988,41(3):281-293.
[9]Veers P S. Modeling stochastic wind loads on vertical-axis turbines[R]. Technical Report SAND83-1909, Sandia National Laboratories,A1buquerque,New Mexico,1984.
[10]Veers P S. Three-dimensional wind simulation[R]. Technical Report SAND88 -0152, Sandia National Laboratories,A1buquerque,New Mexico,1988.
[11]Salzmann D J C. Dynamic response calculation of offshore wind turbine monopile support structures [D]. Delft University of Technology,Netherlands,2004.
[12]Tempel J V D. Design of support structures for offshore wind turbines[D]. Delft University of Technology,Netherlands,2006.
[13]賀廣零,李杰.風力發電機組風致動力響應分析[J]. 電力建設,2011,32(10):1-9.
[14]賀廣零,李杰.基于物理機制的風力發電高塔系統風場模擬[J].同濟大學學報:自然科學報,2010,38(7):976-981.
[15]Li J,Chen J B. Stochastic dynamics of structures [M]. John Wiley &Sons,2009:11-32.
[16]賀廣零,李杰. 風力發電高塔系統基于物理機制的旋轉樣本功率譜研究[J]. 中國電機工程學報,2009,29(26):85-91.
[17]賀廣零,李杰. 風力發電高塔系統風致隨機動力響應分析[J].振動工程學報,2011,26(4):696-703.
[18]賀廣零,李杰. 風力發電高塔系統抗風動力可靠度分析[J]. 同濟大學學報:自然科學報,2011,39(4):474-481.