賀廣零,李 杰

(同濟大學 建筑工程系,上海200092)

隨機風場作用下結構的響應以及動力可靠度是工程界普遍關心的問題.為了分析結構的時域風振響應,需要通過風場模擬,生成具有結構所在場地風荷載統(tǒng)計特性的樣本函數.研究現有的風場模擬方法(如諧波合成[1-2]、線性濾波器[3-4]和小波分析[5-6]等),不難發(fā)現,這些方法都基于功率譜密度函數.而在本質上,功率譜密度函數(包括自譜和互譜)是對平穩(wěn)隨機過程的數值特征描述,很難反映本源隨機過程的豐富概率信息,由此也導致了對結構進行精細化隨機振動反應分析的困難[7].事實上,即使是平穩(wěn)過程,僅依據數值特征解答也很難獲得結構動力可靠度的精確解答.近期,李杰和張琳琳從隨機函數的角度出發(fā)研究隨機過程,將脈動風速時程樣本經Fourier 變換得到的Fourier 譜視為隨機函數,提出了隨機Fourier 譜的概念, 獲得了隨機Fourier 幅值譜經驗物理公式,從而探索了一條從物理本質上反映隨機過程的新道路[8].在隨機Fourier譜的基礎上,賀廣零和李杰針對風力發(fā)電高塔系統(tǒng),提出了考慮槳葉旋轉效應的旋轉Fourier 物理譜[9].在此基礎上,筆者首先提出了基于物理機制的旋轉Fourier 互譜,有效地考慮了槳葉旋轉效應和槳葉上不同點風速之間的相關性.結合典型的1.25 MW 三槳葉變槳距風力發(fā)電高塔系統(tǒng),對旋轉Fourier 譜(隨機Fourier 譜)進行逆Fourier 變換,生成了作用于槳葉(塔體)的風速時程.

1 旋轉Fourier 譜物理機制

空間中一點的風速,可以用一個標準的隨機過程描述.對于體量較大的結構,需要反映不同位置風速之間的相互關系,此時空間的風速是一個時變的隨機場.對于時變隨機場,空間中任意一點的風速不僅與該點的空間位置坐標有關,而且與時間有關.

為了準確描述作用在旋轉槳葉上的風速時程,可有規(guī)律地在風輪平面取樣(圖1).假設在同一半徑r上取N個樣本點,并將風輪平面上l1點風速時程中零時刻的風速……l i點風速時程中(i-1)Δt時刻風速,l i+1點風速時程中iΔ 時刻風速……l N點風速時程中(N-1)Δt時刻的風速提取出來,按照時間順序組合成一組新的風速時程,對其進行Fourier 變換即可得到旋轉Fourier 譜.在本質上,旋轉Fourier 譜由作用在槳葉上的風速經過Fourier 變換而來,是一種自身蘊含了槳葉旋轉效應的紊流風速譜.因此,基于該譜進行風力發(fā)電高塔系統(tǒng)風振動力響應分析時,無須再次考慮槳葉旋轉效應.

圖1 風輪平面上取樣點位置Fig.1 Sample points in the rotor plane

2 旋轉Fourier 互譜模型

槳葉上2 點處脈動風速的旋轉Fourier 互譜,與塔體上2 點處脈動風速的隨機Fourier 互譜有本質的不同.主要體現在兩個方面:①由于槳葉的旋轉效應,旋轉Fourier 互譜必須在旋轉坐標系下考慮2點處脈動風速的相關性.在旋轉坐標系下,2 點處脈動風速的互譜已經不能簡單地通過各點處脈動風速的自譜與相干函數的乘積來確定.②旋轉Fourier互譜體現在槳葉(而非風輪平面)上任意2 點處脈動風速之間的相關性.因此,旋轉Fourier 互譜可分為同一槳葉上2 點之間的旋轉Fourier 互譜和不同槳葉上2 點之間的旋轉Fourier 互譜.

已知槳葉以轉速n0勻速旋轉,半徑r i上的一點在t時刻的風速時程幅值為x i,半徑r j上的一點在t+τ時刻的風速時程幅值為x j(圖2).xi和x j的互相關函數可定義為

x i和x j的Fourier 互譜可表示為[10]

式中,Fxi(n),F xj(n)為隨機Fourier 譜,根據隨機過程的隨機函數描述,可定義為[10]

式中:X(η,t)為隨機過程樣本x(η,t)的集合;η為影響隨機激勵發(fā)展過程且具有物理意義的隨機變量或隨機向量.文獻[11] 由310 組實測風速數據記錄識別出隨機變量地面粗糙度z0服從對數正態(tài)分布,10 m 高平均風速v10服從極值Ⅰ型分布,并最終確定隨機Fourier 譜表達式為

式(2)中,γ(d(τ))為相干函數,其表達式[12]為

式中:a為衰減常數,一般通過實驗獲得,文獻[13]建議取為10;n為頻率;vh為輪轂處平均風速;d(τ)為2 點之間距離(圖2).其表達式為

式中,φ為相位因子.當i,j這2 點處于同一片槳葉時,φ=0;處于不同槳葉時,φ則為2π/ Nb(Nb為槳葉數目)的整數倍.對于三槳葉風力發(fā)電高塔系統(tǒng)而言,φ始終為2π/3.

Fourier 互譜與互相關函數構成Fourier 變換對.對式(2)逆Fourier 變換,可得兩點間的互相關函數為

根據旋轉Fourier 譜的物理機制,旋轉槳葉上任意2 點在不同時刻的旋轉互相關函數,可用2 點之間的互相關函數來代替

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

將式(7),(8)代入式(9)中,

為了分析方便,可將γ(d(τ),n′)進行Fourier 展開

式中:n0為槳葉轉動頻率;k m(n)為Fourier 展開系數,即

將式(11)代入式(10)中,

亦即

(14)可改寫為

進而可得

當r i=r j,即相位因子φ=0 時,旋轉Fourier互譜退化為旋轉Fourier 自譜

式(17)與文獻[9]中的結果是完全一致的.這一結果說明,旋轉Fourier 互譜與旋轉Fourier 自譜具有統(tǒng)一性.前者具有一般性,而后者是特殊情況下的旋轉Fourier 互譜.同時也表明,旋轉Fourier 互譜不再是旋轉Fourier 自譜與相干函數的乘積.

3 隨機風場模擬

3.1 脈動風速

譜和時程這兩類描述方法具有等價性,均可以用來描述隨機過程.只不過前者是對隨機過程的頻域描述,后者是對隨機過程的時域描述.本質上,樣本Fourier 譜和時程共同構成Fourier 變換對.因此,在隨機變量給定后,隨機函數模型轉化為確定性函數,對其逆Fourier 變換即可獲得風速時程.基于這一核心思想,以下完成風力發(fā)電高塔系統(tǒng)風場仿真.

由旋轉Fourier 譜和旋轉Fourier 互譜,可共同構造一維多變量零均值隨機過程的隨機Fourier 譜矩陣如下:

式中,對角線元素由旋轉Fourier 幅值譜組成,非對角元素由旋轉Fourier 互譜組成.需要說明的是,通常情況下互譜總是復數形式的,但是在風工程中,認為在大氣中互譜密度函數的虛部(即正交譜)相比其實部(即互譜)是很小的,工程應用上可以忽略不計[14].從而,這里的旋轉Fourier 互譜也就相應地成為了實數的形式,這是在風工程中出現的特例.

為模擬風場, 首先對隨機Fourier 譜矩陣進行Cholesky 分解

I(n)為下三角矩陣,且有如下形式:

由于隨機Fourier 譜矩陣是n的實值偶函數矩陣,所以上式中元素I ij(n)=I ij(-n).根據式(19),將分解后,就可以用下式對目標隨機過程v j(t),j=1,2, …,k進行仿真,即

式中:φ0ml為隨機初相位角,在[0,2π]區(qū)間取值;Δφml為相位差譜[15];n ml為雙索引頻率,按下式取值:

式(23)中,n u為截斷頻率.由于一般隨機Fourier 譜函數的頻率分布區(qū)間為無窮大,為了數學上處理的方便,有必要設置截斷頻率n u,認為超過該上限后的隨機Fourier 譜函數值為0.由于存在截斷頻率n u, 可知當M→∞時, Δn→0, 因此, 有n u=MΔn.截斷頻率的值由在區(qū)間[0,n u]和區(qū)間[0, ∞]中隨機Fourier 幅值譜下包含的面積之比來確定,通常要求該比例接近于1.可按下述公式表示該準則:

同時,為了防止混疊,根據采樣定理,在使用式(21)生成時程樣本時,時間步長Δt應滿足如下條件

塔體脈動風速生成與槳葉相似,只需將旋轉Fourier 譜換成隨機Fourier 譜即可.限于篇幅,茲不贅述.

3.2 平均風速

對于風力發(fā)電高塔系統(tǒng)而言,風剪模型通常采用指數模型

式中:vh為輪轂高度處的平均風速;zh為輪轂高度;α為風速廓線指數.

對于旋轉槳葉上的任意一點,其高度z因槳葉旋轉而呈現周期性變化

式中:r為計算半徑,指風輪旋轉平面內任意一點與輪轂中心之間的距離;φ=Ωt,為該點在風輪平面的方位角,正上方時為0°, Ω為槳葉旋轉速度.

將式(27)代入式(26)中,可得槳葉上半徑r處的平均風速vs(r)為

4 數值算例

以典型的1.25 M W 三槳葉變槳距風力發(fā)電高塔系統(tǒng)為例,進行風場仿真研究.該風力發(fā)電高塔系統(tǒng)輪轂高度為68 m,風輪直徑為64 m,槳葉轉速為21.1 r ·min-1(0.352 Hz).根據有限元方法,對風力發(fā)電機整體結構進行離散,每片槳葉等效為3 個均勻分布的集中質點,三片槳葉一共為9 個集中質點.塔體(機艙)等效為非均勻分布的6 個集中質點,各點的具體位置依據計算方便原則選定,其簡化動力計算模型如圖3 所示.等效集中質點為動力計算時需要輸入風速時程的計算點,這里主要模擬這些點上的風速時程.

圖3 風力發(fā)電高塔系統(tǒng)計算點Fig.3 Computing points of the wind turbine system

根據式(24),可確定截斷頻率nu.結合實際需要,這里取nu=10 Hz.總頻點數取為M=2 048.繼而由式(23)可得,Δn=0.004 88 Hz.編制Matlab 程序仿真脈動風速時程.表1 給出了塔體第10 ～15 點處的平均風速.圖4 給出了塔體第10 ～15 點的仿真脈動風速時程.不難發(fā)現,不同點的風速時程之間存在一定的相關性,且相關程度隨著2 點距離的增加而減少.例如,相鄰點風速時程之間的相似程度要大于非相鄰點.

表1 塔體各計算點的平均風速Tab.1 Mean wind velocities at the computing points of the wind turbine tower

圖4 塔體各計算點的仿真脈動風速時程Fig.4 Fluctuating wind velocities at the computing points of the wind turbine tower

圖5 給出了槳葉第1 ～3 點處的平均風速時程,圖6 給出了3,6,9 點處的平均風速時程.總體上,旋轉槳葉上各點的平均風速具有如下特點:①平均風速不再為定值,而呈諧波規(guī)律變化;②計算點半徑越大,風速波動幅度越大,如點2 的波動幅度大于點1,點3 的波動幅度大于點2;③不同槳葉之間風速不同步,相鄰槳葉之間存在2π/nb(nb為槳葉數目)的相位差.槳葉1,2,3 要落后于槳葉4,5,6,而4,5,6落后于7,8,9,它們之間的相位差均為2π/3.

圖7 給出了槳葉第1 ～3 點處的脈動風速時程.圖8 給出了固定點風速時程(基于隨機Fourier 譜的風速時程)與旋轉點風速時程(基于旋轉Fourier 譜的風速時程)之間的比較.相比較而言, 基于旋轉Fourier 譜的風速時程具有兩個基本特點:①風速時程幅值有一定增大,但不是特別顯著;②風速時程振動頻率有大幅度提高.這點極為顯著.造成這些差別的主要原因在于旋轉點風速時程不僅體現了風速自身的脈動特性,而且還刻畫了旋轉槳葉高度周期性變化引起的風速波動.總體上,塔體上某點的脈動風速具有時變性,而旋轉槳葉上某點的脈動風速時程具有時變、空間變化雙重特性.在考慮槳葉旋轉效應之后,槳葉脈動風速時程幅值存在一定增長,振動頻率會有大幅度提高,從而必然會對風力發(fā)電高塔系統(tǒng)極值荷載和疲勞荷載產生重要影響.這正是研究槳葉旋轉效應、提出旋轉Fourier 譜的根本意義.

圖5 同一槳葉上不同點的平均風速時程比較Fig.5 Mean wind velocities at different computing points of the same blade

圖6 不同槳葉上半徑相同點的平均風速時程比較Fig.6 Mean wind velocities at different computing points with the same radius

圖7 計算點1, 2,3 處的脈動風速時程Fig.7 Fluctuating wind velocities at computing points 1,2 and 3

圖8 旋轉點風速時程與固定點風速時程比較(100～150 s)Fig.8 Comparison between the fluctuating wind velocities of the rotating point and the stationary point(100～150 s)

圖9 給出了槳葉第3 點處仿真脈動風速時程的計算旋轉Fourier 幅值譜和目標旋轉Fourier 幅值譜的比較圖,其中,旋轉Fourier 幅值譜是由模擬生成的10 000 次、持時為600 s 的風速時程的Fourier 變換幅值經平均計算得到,而目標旋轉Fourier 幅值譜是將T=600 s,z0=0.029 m 和v10=14.71 m ·s-1代入表達式(16)得到.其他各點的分析結果與第3點相似,茲不贅述.結果表明,仿真脈動風速時程符合結構所在場地風荷載統(tǒng)計特性,滿足了模擬仿真的要求.

圖9 點3的計算旋轉Fourier 譜與目標旋轉Fourier 譜比較Fig.9 Comparison between the numerical rotational Fourier spectrum and the target rotational Fourier spectrum(point 3)

5 結語

基于物理機制的旋轉Fourier 互譜,準確反映了作用在旋轉槳葉上風速隨機過程的物理本質,有效地考慮了槳葉的旋轉效應,并體現了槳葉上不同點風速之間的相關性,為準確確定槳葉風荷載奠定了基礎.

實現了風力發(fā)電高塔系統(tǒng)的風場模擬.風力發(fā)電高塔系統(tǒng)風場模擬可分為槳葉、塔體兩部分.其中,槳葉風場模擬需考慮槳葉旋轉效應.依據本文風速譜的建模機制,槳葉(塔體)風場模擬本質上為旋轉Fourier 譜(隨機Fourier 譜)的逆Fourier 變換過程.研究表明,計算旋轉Fourier 幅值譜與目標旋轉Fourier 幅值譜吻合良好,本仿真算法可以合理地模擬作用于風力發(fā)電高塔系統(tǒng)的風速時程.

[ 1] Rice S O.Mathem atical Analysis of Random Noise[ C] ∥Selected Papers on Noise and Stochastic Processes.Dover:[ s.n.] ,1954:133-294.

[ 2] Schueller G I,Shinozuka M.Stochastic methods in structural dynamics[ M] .Dordrecht:Martinus Nijhoff Publishers,1987.

[ 3] Iannuzzi A,Spinelli P.Artificial wind generation and structural response[ J] .Journal of Structural Engineering, ASCE, 1987,113(12):2382.

[ 4] Naganuma T , Deodatis G, Shinozuka M.ARM A m odel for tw o-dim ensional process [ J ] .Journal of Engineering Mechanics,ASME,1987,113(2):234.

[ 5] Yamada M, Ohkitani K.Orthonormal w avelet analysis of turbulence[ J] .Fluid Dy namics Research,1991(8):101.

[ 6] Kitagaw a T,Nomura T .A wavelet-based m ethod to generate artificial w ind fluctuation data[ J] .Journal of w ind engineering and industrial aerodynamics,2002,90:943.

[ 7] 李杰,陳建兵.概率密度演化方程——歷史、進展與應用[ R] .上海:同濟大學土木工程學院建筑工程系,2009.LI Jie, CH EN Jianbing.Probability density evolution equation—history , development and application [ R ] .Shang hai:Tongji University.Department of Building Engineering,2009.

[ 8] 李杰,張琳琳.脈動風速功率譜與隨機Fourier 幅值譜的關系研究[ J] .防災減災工程學報,2004,24(4):363.LI Jie,ZHANG Linlin.A study on the relationship between turbulence pow er spectrum and stochastic Fourier amplitude spectrum [ J] .Journal of Disaster Prevention and Mitigation Engineering,2004,24(4):363.

[ 9] 賀廣零,李杰.風力發(fā)電塔基于物理機制的旋轉Fourier 譜[ R] .上海:同濟大學土木工程學院建筑工程系,2008.H E Guangling, LI Jie.Rotational Fourier spectrum of wind turbine systemsbased on phy sical mechanism [ R] .Shanghai:Tongji University.Department of Building Engineering,2008.

[ 10] 李杰.工程結構隨機動力激勵的物理模型[ R] .上海:同濟大學土木工程學院建筑工程系,2008.

LI Jie.Physical stochastic models for the dy namic excitations of engineering structures [ R] .Shanghai:Tong ji University.Department of Building Engineering,2008.

[ 11] 李杰,張琳琳.實測風場的隨機Fourier 譜研究[ J] .振動工程學報,2007,20(1):66.LI Jie, ZHANG Linlin.Research on the random Fourier spectrum of observational wind [ J] .Journal of Vibration Engineering,2007,20(1):66.

[ 12] World Meteorological Organization.M eteorological aspects of the utilizationof wind as an energy source [ R] .WMO Technical Note No.175 World M eteo rological Organization,Geneva 1981.

[ 13] Drag t J B.T he spectra of w ind speed fluctuations m et by a rotating bladeand resulting load fluctuations[ C] ∥Proceedings of European Wind Energy Conf 1984.H am burg:[ s.n] , 1984:453-459.

[ 14] Simiu E, Scanlan R H .Wind effects on structures:an introduction to w ind engineering [ M ] .New York:John Wiley,1978

[ 15] 艾曉秋.基于隨機地震動模型的地下管線地震反應及抗震可靠度研究[ D] .上海:同濟大學建筑工程系,2005.AI Xiaoqiu.Stochastic earthquake model-based research on seismic responseand reliability of underground pipelines [ D] .Shang hai:Tongji University.Department of Building Engineering,2005.