脈動風速時程數值模擬

程　華，鐘華生，周　凌，王　龍

(后勤工程學院　a.軍事土木工程系; b.學員旅，重慶　401331)

?

脈動風速時程數值模擬

程華a，鐘華生a，周凌a，王龍b

(后勤工程學院a.軍事土木工程系; b.學員旅，重慶401331)

摘要：對高層、高聳和大跨輕柔等風敏感結構而言，風荷載是主要控制荷載，必須考慮風荷載的動力作用。分析了風速時程模擬方法的概況，包括傳統的諧波疊加法、線性濾波法，以及近年來發展起來的小波分析、逆傅立葉變換法。采用AR法模擬了某索-桁移動機庫的脈動風速時程?？偨Y了各種方法的優缺點，展望了風速時程模擬的前景。

關鍵詞：脈動風速時程；數值模擬；諧波疊加法；線性濾波法；小波分析

隨著工程結構日趨多樣化、大型化、復雜化，新材料、新工藝、新技術、新形式以及新設計方法的應用，使得結構更加輕柔，結構體型更加復雜，從而使得結構對于風荷載的作用更加敏感，風荷載已成為結構抗風計算、防災減災分析的重要設計荷載，結構風振響應研究日益受到重視[1]。目前，雖然在一些較為重要建筑中采用風洞試驗或現場實測等手段確定結構承受的風荷載，但這些方法都比較復雜，耗時耗資，不具備通用性。此外，我國地域廣闊，地區之間的風速特性也很不相同，已有的強風記錄資料有限，不能在結構風振分析中普遍應用。因此，風荷載數值模擬方法在實際工程中得到了廣泛采用。

1諧波疊加法

諧波疊加法是基于一系列三角級數求和，采用以離散譜逼近目標隨機過程的一種平穩隨機過程數值模擬方法。在1954年Rice[2]首先提出諧波疊加法的基本思想，但當時僅適用一維單變量平穩高斯隨機過程的模擬。在1972年Shinozuka[3]提出CAWS法和WAWS法，將諧波疊加法發展到模擬多變量、非平穩高斯隨機過程。雖然后來Yang等[4]人在諧波疊加法中引入快速傅立葉(FFT)變化技術，大大加速了諧波疊加法的模擬速度，但它的計算量依然較大。采用諧波疊加法模擬風荷載時，需要在每個頻率上進行大量運算，并且不能考慮時間相關性。

隨機信號可以通過離散Fourier變換分析，分解為一系列具有不同頻譜和幅值特性的正弦或其他諧波，譜密度就等于由帶寬劃分的這些諧波幅值的平方。根據Shinozuka的理論[3]，一組m個s維隨機過程xj(t)(j=1,2,…,m)可模擬為：

(1)

風速時程的互功率譜密度函數矩陣一般形式為：

(2)

(3)

對Sij(n)進行喬里斯基(Cholesky)分解：

(4)

可以得到下三角矩陣元素H(n)，H*(n)Τ是H(n)的轉置復共軛矩陣。

θjk(nl)表示Hjk(nl)的相位角：

(5)

為避免模擬結果失真，時間間隔應滿足Δt≤π/(2nl)。

當模擬的點大于200時，按常規方法計算最終的風速時程是很耗機時的，建議引入FFT技術[6]，將大大提高計算效率。

取M=2π/(Δt·Δn)為整數，則可將式(1)改寫為：

(6)

Gj(p·t)可用傅立葉變換計算：

(7)

其中：

(8)

諧波疊加法簡單直觀、數學基礎嚴密，很多學者采用經典諧波疊加理論或結合快速傅立葉變換技術模擬脈動風荷載[7-11]，做了大量工作，但很少取得新的突破。李正農等[12]采用非均勻圓頻率間隔的方法，改進了以往諧波合成法模擬的周期性缺點，但此種方法是以犧牲效率為代價的。Debasis等[13]在分析總結高斯模型的基礎上，提出了一種新的非高斯條件的迭代方法，結合非線性非高斯映射技術，模擬了大跨橋梁的風速時程曲線；并用不同模擬方案的風荷載對大橋進行風振響應分析，與實測數據對比表明，非高斯模擬方案優于高斯模擬方案。Huang等[14]提出了一種全新的針對功率譜密度函數喬里斯基(Cholesky)分解方法，它的關鍵思想是從復雜EPSD/PSD矩陣分解中分離出各個獨立的階段，最后再以某種特定的方式進行組裝，進而大大提高分解效率。

2線性濾波法

線性濾波法則是基于數字濾波技術，將均值為零的白噪聲隨機序列通過線性濾波器，輸出具有指定譜特征的隨機過程。線性濾波器包含AR法、MA法和ARMA法，其中AR法被廣泛應用于脈動風速時程模擬。

脈動風速時程本質上也是隨機時間序列，因此可以采用AR法模擬脈動風的風速時程。m個空間相關脈動風速時程v(X,Y,Z,t)列向量的AR模型可表示為：

(9)

其中：X=[x1,x2,…,xm]Τ、Y=[y1,y2,…,ym]Τ、Z=[z1,z2,…,zm]Τ；(xi,yi,zi)為空間第i點的坐標，i=1,2,…,m；Δt為模擬風速時程的時間步長；p為模型的階數；φk為AR模型的自回歸系數矩陣，為m×m階矩陣，k=1,2,…,p；N(t)為獨立隨機過程向量。

為簡便起見，以下將v(X,Y,Z,t)簡記為v(t)，將式(9)的左右兩邊同乘vΤ(X,Y,Z,t-jΔt)，并同時取數學期望可得相關函數R(jΔt)與回歸系數φk的關系為：

(10)

將式(10)改為矩陣形式為：

(11)

其中：[R]pm×m=[R(Δt),R(2Δt),…,R(pΔt)]Τ，[φ]pm×m=[φ1,φ2,…,φp]Τ。

根據隨機振動理論，相關函數可由維納-辛欽公式求出

(12)

式中：n為脈動風速頻率；Sij(n)為風速譜函數，可由脈動風速自譜密度函數Sii(n)和相干函數cohij(n)(i=1,2,…,m; j=1,2,…,m)確定

(13)

獨立隨機過程向量N(t)由下式求解

(14)

式中：n(t)=[n1(t),n2(t),…,nm(t)]Τ，ni(t)是零均值、方差為1且彼此相互獨立的正態隨機過程，i=1,2,…,m；L為m階下三角矩陣，通過m×m階協方差矩陣RN的喬里斯基(Cholesky)分解確定

(15)

計算時，假定初始時刻之前的風速為0，即： t≤0時,v(t)=0。最終的m個隨機脈動風速時程可表示為

(16)

文獻[15-17]中采用基于AR模型的線性濾波法對超高層結構的脈動風速時程進行了數值模擬，驗證了其可行性與有效性，且占用內存小、計算速度快。馬駿等[18]結合諧波疊加法和自回歸法的優點，提出風速時程的高效高精度新型混合法。在保證計算精度下，該方法的模擬效率是諧波疊加法的數十倍；計算精度與自回歸方法相比亦顯著提高，且易于收斂。曾少青對自回歸階數P=4，5，6，7時的風速譜模擬效果進行對比，發現P=6時模擬譜與目標譜的標準差最小，說明并非階數越高風速譜模擬精度越高[19]。武秀根等采用AR模型模擬了超限冷卻塔的表面風壓表明，基于隨機脈動風荷載模擬的數值計算方法，特別適用于冷卻塔高度超出規范要求的情況下評估冷卻塔的風振響應[20]。劉亞慶分析了AR法結合快速傅立葉變換的頻率不準確采樣誤差，由此提出了兩種改進方法：一是精確頻率小數位數，自振頻率的計算結果每多保留一位就會使采樣點的個數變為原來的10倍，計算量巨大；二是構造結構相鄰兩階自振頻率的差值向量，使其模最小，取得了較為顯著的效果[21]。

3小波分析

設φ(t)為一平方可積函數，即φ(t)∈L2(R)，L2(R)為平方可積的函數空間，如果其傅立葉變換φ(ω)滿足條件[22]

(17)

則稱φ(t)為基本小波或小波母函數，而連續小波基函數φα,τ(t) 則通過小波母函數φ(t)的伸縮和平移得到。由同一母函數φ(t)經過伸縮和平移后得到一組函數系列

(18)

式中：α為尺度因子；τ為平移因子。

風工程中使用的信號大都為離散函數，因此對小波加以離散化。在數值計算中最常用的是二進制采樣網格。每個網格點對應的尺度為2j，平移量為2jk的小波稱為二進小波(Dyadic Wavelet)[23]，即：

(19)

(20)

相應可以將隨機信號X(t)表示為一系列離散小波基的疊加[24]

(21)

其中，小波系數Wj,k與尺度系數Uj,k可表示為[25]

(22)

由于自身良好的時頻特性，小波分析受到了國內外學者的廣泛關注。吳筑海等[26]用分解后的小波系數來描述風速信號在時間與頻率兩個域上的特性，證實了小波分析在風速時程分析中的穩定性。韓艷等[27]采用正交Mayer小波變換模擬脈動風，發現模擬的脈動風具有一定的間歇性和局部相似性，Yamada等[28]對實測的脈動風進行小波分析，得到了類似結論。Zeldin等[29]采用db3小波，利用AR法來合成每一尺度上的小波系數，而采用多個小波函數來定義多維隨機風場模擬。Aksoy等[30]采用Haar小波模擬了1000年風速，并與5種傳統的參數模型模擬方法分別進行了比較，發現小波是保留信號相關性的最好方法。

陳艾榮等[31]根據小波的多分辨率分析性質，利用小波重構算法，選取db3小波基，基于線性估計原理，利用逆小波變換反復迭代計算，模擬了蘇通橋某點的縱向脈動風速。周岱等[32]采用離散小波變換分解和重構風速時程，并引入最大熵值法減少風速時程分析在時-頻域上的信息損失并提高效率。柳杰等[33]采用Haar小波對風速時程進行分解和識別，提高分解后的低頻分量，使得風速時程在時-頻域上更加逼近目標值。郭立偉等[34]利用db3小波，重構14層，根據AR模型模擬出21條具有空間相關性的脈動風速時程。

李宏男等[35]系統總結了小波分析在結構風工程中的風荷載模擬、風場特性研究、風與風壓的相互關系以及模態參數識別中的國內外現狀，并對小波變換在結構風工程中的應用前景進行了展望。

4逆傅立葉變換法(IFFT法)

逆傅立葉變換法(IFFT法)[36-37]由D.Cebon[38]首先提出，在國內最早用于鐵路不平順隨機過程的模擬[39]。候傳亮等[40]針對道路過程，引入道路車輛速度后轉換空間頻率為時間頻率，運用IFFT法和諧波疊加法成功模擬了一種C級道路時間序列。逆傅立葉方法(IFFT)的基本思想是通過時間序列估計功率譜密度的周期圖法—Blackman-Tukey法來反推離散后的功率譜密度(PSD)與時間序列的關系：

(23)

IFFT法為風速時程的數值模擬提供了一種新思路。根據已知的風速功率譜密度，則可以由式(24)反算隨機過程的傅立葉頻譜值。進而通過對復序列F(k)的IFFT的時序化模型

(24)

5算例分析

以某30 m×40 m索-桁機庫(圖1)為工程背景，該機庫主體結構由10榀索-桁架(圖2)組成，兩端設有蚌殼式門，骨架采用輕型鋁合金型材，屋面敷設M950-2型膜材。采用基于AR模型的線性濾波法，運用Matlab編程模擬脈動風時程，模擬的具體參數見表1所示。

圖1　機庫結構簡圖

圖2　一榀索-桁架簡圖

圖3表示生成的第1點和第10點的脈動風速時程曲線；圖4采用雙對數坐標，給出了模擬脈動風速自譜功率譜密度與目標譜的對比曲線?？梢钥闯?，脈動風速在0上下波動，符合平穩隨機過程的特性；模擬風速譜和目標譜吻合較好(尤其在高頻階段)，表明數值模擬的脈動風速時程自相關隨機特性滿足所需隨機特征。

表1　風速時程模擬時的具體參數

圖3　模擬脈動風速時程

圖4　模擬脈動風速功率譜分析

6結束語

本文介紹了模擬脈動風速時程的幾種方法，重點介紹了諧波疊加法和基于AR模型的線性濾波法的基本原理，討論了各種方法的研究近況，并采用AR法模擬了脈動風速時程。

諧波疊加法和線性濾波法都是模擬平穩高斯隨機過程的經典方法，能夠快速高效地滿足工程精度要求，尤其是引入快速傅立葉變換技術后，大大提高了計算效率，在工程界和學術界得到了廣泛應用。但諧波疊加法要求諧波頻率必須滿足均勻分布的條件，且僅適用于線性情況，也不能考慮風速的時間相關性；線性濾波法算法繁瑣、精度較差，AR模型定階還值得深入研究。

小波分析被譽為“數學顯微鏡”，在時域和頻域中有良好的局部化特征，它能夠聚焦到風速時程的任意細節并加以分析，快速準確地提取樣本的局部譜密度特征，并可用局部能量密度函數表示風頻率隨時間的變化。作為一種新的模擬工具，小波分析在風速時程的模擬研究中具有廣闊發展前景。

目前，風速的模擬與結構形式毫不相關，也就是說不論什么樣的結構形式，只要所在的場地地貌相同，平均風速相同，模擬出來的脈動風速時程都是相同的，這與實際情況存在一定差異。脈動風速與結構形式的關聯值得進一步研究。

參考文獻：

[1]MARIO D P.Digital Simulation of Wind Field Velocity[J].Journal of Wind Engineering and Industrial Aerodynamics,1998,74-76:91-100.

[2]RICE S O.Mathematical analysis of random noise[C]// Selected papers on noise and Stochastic processes,N.Wax,et.Dover Publish Inc.New York,1954:133-294.

[3]SHINOZUKA M.Simulation of multivariate and multimensional random processes[J].Journal of the Acoustical Society of America.1971,49(1):357-367.

[4]YANG J N.Simulation of random envelope processes[J].Journal of Sound and Vibration,1972,21(1):73-85.

[5]王之宏.風荷載的模擬研究[J].建筑結構學報，1994，l5(1)：44-52.

[6]司敏，陳波，劉元志，等.斜拉人行橋多維脈動風荷載模擬研究[J].武漢理工大學學報，2014，36(8)：82-87.

[7]李志國，何能，廖海黎，等.數值模擬脈動風場在斜拉橋抖振響應分析中的應用[J].鐵道建筑，2014(5)：9-12.

[8] 汪瑞.結合工程實例的扣塔脈動風荷載模擬及風振響應分析[J].特種結構，2012，29(4)：124-127.

[9]曹玉生，包格瓦日圖.基于Matlab的大型兆瓦級風電機脈動風速時程數值模擬[J].內蒙古工業大學學報，2013，32(4)：278-284.

[10]舒興平，嚴興華，鄒浩.通信角鋼塔脈動風荷載模擬及風振響應分析[J].工業建筑，2011，41(4)：119-123.

[11]曹麗珍.高聳鋼塔結構風振響應分析[J].四川建材，2014，40(4)：59-62.

[12]李正農，劉艷萍，王鶯歌，等.近地邊界層脈動風數值模擬[J].湖南大學學報:自然科學版，2010，37(4)：6-10.

[13]DEBASIS K,SAMIT R C,MASANOBU S.Probabilistic Engineering Mechanics[J].2012(29):167-175.

[14]HUANG G Q,LIAO H L,LI M S.New formulation of Cholesky decomposition and applications in stochastic simulation[J].Probabilistic Engineering Mechanics,2013(34):40-47.

[15]陳俊儒，呂西林.上海中心大廈脈動風荷載模擬研究[J].力學季刊，2010，31(1)：92-100.

[16]胡亞琪.輸電塔架結構風振響應分析[D].廣西：廣西大學，2011.

[17]陳建芳，陳家豪，何志.基于AR法的輸電塔風振響應時程分析[J].廣西工學院學報，2013，24(2)：74-77.

[18]馬駿，周岱，李磊，等.風時程模擬的高效高精度混合法[J].工程力學，2009，26(2)：53-59.

[19]曾少青.寧波國際會展中心屋蓋結構風振響應分析[D].廣州：華南理工大學，2010.

[20]武秀根，鄭百林，劉曙光，等.大型冷卻塔的風致振動數值分析[J].力學季刊，2011，32(1)：91-96.

[21]劉亞慶.Kiewitt型索彎頂結構的靜力性能與風振響應分析[D].廣州：華南理工大學，2013.

[22]CHARLES K C.An Introduction to Wavelets[Z].Academic Press Inc.1992.

[23]ASSMOAH F.Discrete Wavelet Analysis of Signal[J].International Journal of Engineering Education，1999 (36):255-264.

[24]李峰，竇鵬，王磊.風荷載模擬方法探討[J].水利與建筑工程學報，2009，7(4)：57-59.

[25]KAREEM A,KIJEWSKI T.Time-frequency analysis of wind effects on structures[J].Journal of wind engineering and industrial aerodynamics，2002(90):1435-1452.

[26]吳筑海，馬駿，周岱，等.空間結構風時程模擬及分析[C]//第四屆全國現代結構工程學術研討會論文集.寧波：[出版地不詳]，2004.

[27]韓艷，陳政清.利用小波逆變換模擬隨機風場的脈動風[J].振動工程學報，2007，20(1):52-56.

[28]YAMADA M,OHKITANI K.Orthonormal wavelet analysis of turbulence[J].Fluid Dyn.Res.，1991(8):101-115.

[29]ZELDIN B A,SPANOS P D.Random field representation and synthesis using wavelet bases[J].Journal of Applied Mechanics,1996,63(12):946-952.

[30]AKSOY H,TOPRAK Z F,AYTEK A,et al.Stochastic generation of hourly mean wind velocity data[J].Renewable energy,2004(29):2111-2131.

[31]陳艾榮，王毅.基于小波方法的隨機脈動風模擬[J].同濟大學學報:自然科學版，2005，33(4):427-431.

[32]周岱，馬駿，吳筑海，等.空間結構三維風速時程模擬及其小波分析[J].工程力學，2006，23(3):88-92.

[33]柳杰，周岱，黃劍偉，等.大跨空間結構風模擬的改進方法及小波識別[J].振動與沖擊，2006，25(4):21-24.

[34]郭立偉，張文福，孫曉剛，等.電視塔抗風控制[J].東北林業大學學報，2007，35(8):60-62.

[35]李宏男，伊廷華，王國新.小波分析在結構風工程中的應用[J].振動與沖擊，2008,27(9)：11-16.

[36]劉錫良，周穎.風荷載的幾種模擬方法[J].工業建筑，2005，35(5)：81-85.

[37]鄧文旭.結構的脈動風荷載模擬方法探討[J].安徽建筑，2007(4)：85-87.

[38]CEBON D.Interaction Between Heavy Vehicles and Roads[Z].SAE,Inc.1993.

[39]陳果，翟婉明.鐵路軌道不平順隨機過程的數值模擬[J].西南交通大學學報，1999，34 (2):138-142.

[40]侯傳亮，張永林.工程平穩隨機過程的數值模擬研究[J].武漢工業學院學報，2003,22(3)：27-29.

(責任編輯周江川)

本文引用格式：程華，鐘華生，周凌，等.脈動風速時程數值模擬[J].兵器裝備工程學報，2016(4):143-148.

Citation format:CHENG Hua, ZHONG Hua-sheng, ZHOU Ling, et al.Numerical Simulation of Fluctuating Wind Velocity Time Series[J].Journal of Ordnance Equipment Engineering,2016(4):143-148.

Numerical Simulation of Fluctuating Wind Velocity Time Series

CHENG Huaa, ZHONG Hua-shenga, ZHOU Linga, WANG Longb

(a.Department of Civil Engineering; b.Brigade of Cadet,Logistic Engineering University, Chongqing 401331, China)

Abstract:Wind load was the key control load and its dynamic effect for high-rise and long-span soft structures which were wind-sensitive must be considered. The general situation in numerical simulation of fluctuating wind velocity time history was introduced including traditional harmony superposition method, linear filtering method, as well as wavelet analysis, inverse Fourier transform. Then, fluctuating wind velocity time series of cable-truss mobile hangar was simulated in AR model. Finally, the advantages and disadvantages of various methods were summarized and the prospect of wind velocity time history simulation was looked forward.

Key words：fluctuating wind velocity time series；numerical simulation；harmony superposition method; linear filtering method; wavelet analysis

文章編號：1006-0707(2016)04-0143-06

中圖分類號：TU312.1

文獻標識碼：A

doi:10.11809/scbgxb2016.04.035

作者簡介：程華(1958—)，男，工學博士，教授，主要從事軍事工程搶修搶建技術與裝備研究。

基金項目：重慶市基礎與前沿研究計劃項目(cstc2013jcyjA30013)

收稿日期：2015-09-14；修回日期：2015-10-10

【基礎理論與應用研究】