索膜結(jié)構(gòu)風(fēng)致流固耦合氣動阻尼效應(yīng)研究

盧 旦，樓文娟，楊 毅

(1.華東建筑設(shè)計院有限公司，上海 200041;2.浙江大學(xué)土木工程學(xué)院，杭州 310027;3.中國中鐵二院工程集團(tuán)有限責(zé)任公司，成都 610017)

當(dāng)前，建筑物的風(fēng)振計算可分為兩類:即建筑物的自身振動對脈動風(fēng)力變化(產(chǎn)生附加氣動力)有顯著影響和沒有影響。后者可以通過先求得靜止建筑物的風(fēng)力荷載，再進(jìn)行結(jié)構(gòu)的動力響應(yīng)分析的方法進(jìn)行［1］。對于索膜結(jié)構(gòu)等柔性結(jié)構(gòu)而言，由于結(jié)構(gòu)在風(fēng)力作用下的變形相對結(jié)構(gòu)本身尺寸相當(dāng)大時，附加氣動力的作用不可忽略，因此必須采用考慮流固耦合的計算與分析。

流固耦合計算又分為兩種:一種稱之為緊耦合(或強(qiáng)耦合)方法，是對流體與固體的聯(lián)立方程式進(jìn)行一起求解。例如，流體與固體通過有限元方法進(jìn)行離散化后再組合形成一個矩陣進(jìn)行求解。這個方法的主要缺點是不能使用已有的得到廣泛應(yīng)用的流體和結(jié)構(gòu)求解器，而且關(guān)于積分的推導(dǎo)也是非常冗長乏味的。因此，在實際應(yīng)用中，復(fù)雜的流場跟彈性變形較大的結(jié)構(gòu)之間的相互作用的模擬就變得非常費時費力。另一種稱之為松耦合(或弱耦合)方法，即在每一個時間步內(nèi)，通過流體計算得到風(fēng)壓、風(fēng)力并將其輸入到結(jié)構(gòu)動力計算中，求得結(jié)構(gòu)的位移、速度及加速度等動力響應(yīng)。然后，基于這些響應(yīng)在流體計算中對建筑周圍的計算網(wǎng)格進(jìn)行移動從而改變建筑表面的邊界條件，重新計算結(jié)構(gòu)表面的風(fēng)荷載。松耦合方法最大的優(yōu)勢是可以利用現(xiàn)有的通用流體和結(jié)構(gòu)軟件，并且可以分別對每一個軟件單獨地制定合適的求解方法，而流體以及固體的求解器不需要改寫，可以通過聯(lián)合現(xiàn)有的軟件的方法從而得到耦合問題的一個精確的求解。此外，這種方法還可以針對一個CFD程序選擇多個CSD程序。因此可以根據(jù)工程上的不同需求，進(jìn)行流體和結(jié)構(gòu)模型的相互交換，使這種方法更具有一般性。但是，松耦合方法也存在計算效率和收斂性等問題的限制［2］。近年來，隨著功能強(qiáng)大、具有二次開發(fā)功能的商用CFD(Computational Fluid Dynamics)和FEM(Finite Element Mode)軟件以及功能更強(qiáng)大的網(wǎng)格剖分工具的出現(xiàn)，基于松耦合方法的流固耦合數(shù)值模擬技術(shù)的實現(xiàn)逐步成為可能［3-4］。

本文根據(jù)三角形重心坐標(biāo)系原理提出了基于松耦合計算的流固耦合網(wǎng)格協(xié)調(diào)插值新方法，提高了流固耦合的計算效率;同時，為解決由于柔性結(jié)構(gòu)加速度過大而導(dǎo)致的求解不穩(wěn)定問題，通過在1個時間步上進(jìn)行流體及結(jié)構(gòu)間反復(fù)交替求解來得到穩(wěn)定的時程解。最后，通過與風(fēng)洞試驗的對比驗證了本文方法的可靠性。

1 變形跟蹤系統(tǒng)

影響流固耦合計算精度的因素主要有以下幾個方面:流體計算的精度、結(jié)構(gòu)計算的精度和流固耦合界面的插值精度。而流固耦合界面插值精度又包括了節(jié)點力的插值和節(jié)點位移的插值兩大問題。前者通常可以采用“鄰近搜索”等方法實現(xiàn)，而后者往往由于不同網(wǎng)格體系信息量的不對稱導(dǎo)致插值計算的困難［5］。

目前常用的搜索算法主要有以下幾種:鄰近搜索，該方法的困難在于耦合域可能分布于不同程序的好幾個位置上，并且不同的程序的耦合面之間還可能存在著初始距離。強(qiáng)力搜索，即在所有的單元上執(zhí)行一次循環(huán)，顯而易見，這種方法的效率是比較低的。此外，還有八叉樹搜索、超前搜索、桶式搜索等。所有這些算法都基于不同的假設(shè)，并且使用了不同的數(shù)據(jù)結(jié)構(gòu)，并且都具有相應(yīng)的復(fù)雜性，因此給流固耦合的工程應(yīng)用帶來的效率上的困難。

到目前為止，用于流固耦合網(wǎng)格更新的方法很多，如彈簧近似法［6-8］等。這些算法雖然可以保證網(wǎng)格單元的拓?fù)潢P(guān)系不發(fā)生改變，但在處理大變形問題時容易發(fā)生單元交錯，從而導(dǎo)致運算失敗。雖然網(wǎng)格實時重劃［9］計算可以解決上述問題，但該方法耗費機(jī)時，且在網(wǎng)格重畫過程中將破壞原有網(wǎng)格的拓?fù)浣Y(jié)構(gòu)，從而帶來插值誤差。

在索膜結(jié)構(gòu)流固耦合計算中通常使用三角形單元。雖然這些三角形都是空間的，但三角形卻是一個天生的2D物體，原因是三角形的三個點在同一平面上。圖1中三角形代表FEM模型中的一個結(jié)構(gòu)單元，A、B、C分別為結(jié)構(gòu)單元的節(jié)點，點P代表落在三角形ABC內(nèi)部的CFD網(wǎng)格的節(jié)點。假設(shè)以點A作為起點，那么點B相當(dāng)于在AB方向移動一段距離得到，而點C相當(dāng)于在AC方向移動一段距離得到。則對于平面內(nèi)任意一點，可由如下方程來表示:

整理方程(1)得到

令v0=C-A，v1=B-A，v2=P-A則，

等式(3)兩邊分別點乘v0和v1得到兩個等式:

求解上述方程(4)、(5)得到:

如果系數(shù)u或v為負(fù)值，那么相當(dāng)于朝相反方向移動，即BA或CA。因此，如果:u≥0，v≥0，且u+v≤1 則可以判斷點P位于三角形ABC內(nèi)部。上述公式(6)、(7)中僅包含少量的四則運算，因此在計算效率非常高。

基于上述坐標(biāo)系統(tǒng)，原先位于一個結(jié)構(gòu)三角形單元內(nèi)部的流體網(wǎng)格會隨著結(jié)構(gòu)的變形依然貼附于原三角形內(nèi)。并且隨著三角形在平面內(nèi)的變形，三角形內(nèi)部的點僅發(fā)生相對位置的變形，而變形前后的比例關(guān)系保持不變，如圖2所示。

在進(jìn)行流固耦合計算時，流體網(wǎng)格一般總要比結(jié)構(gòu)網(wǎng)格密，即一個結(jié)構(gòu)單元(以三角形單元為例)內(nèi)總會包含多個流體網(wǎng)格節(jié)點。因此，以一個三角形結(jié)構(gòu)單元建立的三角形重心坐標(biāo)系去定位其中所包含的各個流體網(wǎng)格節(jié)點，既能保證各節(jié)點在變形前后的相對位置保持不變，同時又能提高計算效率。

圖1 三角形重心坐標(biāo)系Fig.1 Triangle barycentric coordinate system

圖2 變形跟蹤系統(tǒng)Fig.2 Deformation tracking system

2 收斂性能改進(jìn)

工程界已經(jīng)可以通過使用相對成熟的有限元軟件對結(jié)構(gòu)在風(fēng)荷載下的應(yīng)力、變形等效應(yīng)進(jìn)行精確計算。同時，也能夠通過使用計算流體動力學(xué)軟件對結(jié)構(gòu)所受的風(fēng)荷載進(jìn)行準(zhǔn)確估計。然而，對于進(jìn)行流固耦合計算的方法和工程應(yīng)用卻少之又少。傳統(tǒng)的松耦合在一個時間步內(nèi)僅需進(jìn)行一個循環(huán)的流體計算及結(jié)構(gòu)動力響應(yīng)計算，這種方法對于像空氣與建筑物這樣相對質(zhì)量比大的情況有效。但是，當(dāng)由建筑物變形產(chǎn)生的附加空氣質(zhì)量不能忽略時，以上方法會由于在結(jié)構(gòu)上生成過大加速度從而導(dǎo)致求解不穩(wěn)定。文獻(xiàn)［2］指出可通過在1個時間步上進(jìn)行流體及結(jié)構(gòu)間反復(fù)交替、不斷修正求解，解決以往柔性結(jié)構(gòu)流固耦合計算中的收斂性問題，但未對具體方法進(jìn)行說明。本文以此為思想，提出了一種解決流固耦合計算收斂性的方法，具體過程如下:

圖3所示為傳統(tǒng)松耦合算法。圖中Fi為i時刻結(jié)構(gòu)所受荷載，F(xiàn)i+1為下一時刻i+1時刻需要加載的量。按照傳統(tǒng)算法，即使荷載步被細(xì)分，如圖中所示分為n個荷載步，每次增加dF/n的量，最終和一次加載的效果是一樣的。圖4為改進(jìn)方法，即荷載步被細(xì)分，但Fi+1也隨時進(jìn)行更新，從第二個荷載子步開始，每次增量dF/n變成dF'/n、 dF″/n、dF?/n…以此類推。因此，最終的加載效果如圖4中曲線所示。事實上，這種算法更符合流固耦合的實際發(fā)生機(jī)理，同時也使數(shù)值計算的運行性能得以改善，使計算過程更加容易收斂。

圖3 常規(guī)方法Fig.3 Standard method

圖4 荷載步修正算法Fig.4 Loads step amend arithmetic

3 信號處理方法

3.1 EMD法與Hilbert變換

Huang［10］提出了一種新的信號處理方法—經(jīng)驗?zāi)B(tài)分解方法(Empirical Mode Decomposition，EMD)。Huang［11］又將該方法進(jìn)行了一些改進(jìn)。該方法從本質(zhì)上講是對一個信號進(jìn)行平穩(wěn)化處理，其結(jié)果是將信號中不同尺度的波動或趨勢逐級分解開來，產(chǎn)生一系列具有不同特征尺度的數(shù)據(jù)序列，每一個序列稱為一個本征模函數(shù)(Intrinsic Mode Function，IMF)分量。最低頻的IMF分量通常情況下代表原始信號的趨勢或均值。作為一種應(yīng)用，EMD分解方法可以有效地提取一個數(shù)據(jù)序列的趨勢或去掉該數(shù)據(jù)序列的均值。

EMD分解方法的基本思想是:加入一個原始數(shù)據(jù)序列X(t)的極大值或極小值數(shù)目比上跨零點(或下跨零點)的數(shù)目多2個(或2個以上)，則該數(shù)據(jù)序列就需要進(jìn)行平穩(wěn)化處理。最后，原始的數(shù)據(jù)序列即可由一系列IMF分量以及一個均值或趨勢項表示:

由于每一個IMF分量是代表一組特征尺度的數(shù)據(jù)序列，該過程實際上將原始數(shù)據(jù)序列分解為各種不同特征波動的疊加。

EMD分解的主要目的之一是進(jìn)行Hilbert變換，進(jìn)而得到Hilbert譜。在對每一個IMF分量Cj(t)作Hilbert變換之后，得到一個變換平面內(nèi)的數(shù)據(jù)序列(t):

其中P為Cauchy主值，由Cj(t)和C～j(t)可以構(gòu)成一個復(fù)序列Zj(t):

其中:

得到的瞬時頻率為:

因此，原始數(shù)據(jù)序列可以表示為

3.2 隨機(jī)減量法(RDT)

如果是在環(huán)境激勵下而不是沖擊荷載下采集到隨機(jī)相應(yīng)時程數(shù)據(jù)，在對EMD法所得到的結(jié)果運用Hilbert變換前，必須運用隨機(jī)減量法(Random Decrement Technique，RDT)對各個IMF分量進(jìn)行自由衰減響應(yīng)的提取。該方法通過時間平均，從隨機(jī)振動信號中提取自由衰減響應(yīng)。隨機(jī)減量技術(shù)在使用上簡單易行，物理意義明確，對激勵只有定性要求，且可利用隨機(jī)擾動來激勵，同時只需對響應(yīng)信號進(jìn)行處理，因此該方法在結(jié)構(gòu)模態(tài)參數(shù)識別等許多領(lǐng)域中都得到了成功的應(yīng)用［12-13］。

4 算例分析

本文采用流固耦合松耦合計算方法，聯(lián)合應(yīng)用經(jīng)驗?zāi)B(tài)分解法、隨機(jī)減量法和希爾伯特變換法，對文獻(xiàn)［14］中大跨度柔性屋蓋結(jié)構(gòu)在四周封閉和突然開口兩種情況下的氣動阻尼和氣承剛度進(jìn)行研究，同時探討數(shù)值計算的準(zhǔn)確性。

4.1 風(fēng)洞試驗

本次試驗所采用的是能考慮結(jié)構(gòu)和來流之間相互耦合作用的氣動彈性模型。模型屋面為長寬均等于600 mm的正方形，底裙高200 mm，高跨比H/L=1/3，如圖5所示。采用1 mm厚的鋁合金板作為屋面的材料，彈性模量E=7.1 ×1010N/m2，密度 ρ=2 420 kg/m3，泊松比γ=0.31。同時，為了研究大跨度平屋面在四周封閉和迎風(fēng)面開孔情況下的風(fēng)振響應(yīng)特征，在模型迎風(fēng)面正中開一100 mm×100 mm方孔。對于四周封閉模型，方孔由一塊同等大小的木板封閉。

風(fēng)洞試驗在南京航空航天大學(xué)603研究所的NH-2低速風(fēng)洞中進(jìn)行的。模擬B類地貌大氣邊界層，地貌粗糙度系數(shù)α=0.16。試驗共采用6只加速度傳感器，傳感器在屋面板上的布置如圖6所示。信號的采集與分析由靖江東華測試技術(shù)開發(fā)有限公司研制的DH5935/DH5936動態(tài)信號測試系統(tǒng)完成。封閉模型的采樣頻率為500 Hz，開孔模型的采樣頻率1 000 Hz，采樣持續(xù)時間一般都超過90 s。試驗風(fēng)速為:8 m/s、9 m/s、10 m/s、11 m/s、12 m/s。

圖5 風(fēng)洞模型Fig.5 Wind tunnel model

圖6 加速度傳感器在屋面的測點布置圖Fig.6 Station plan of acceleration sensor on building roof

4.2 流固耦合模型

4.2.1 CFD 模型

在數(shù)值模擬中按1∶1尺寸建立與風(fēng)洞試驗?zāi)Ｐ鸵恢碌挠嬎隳Ｐ汀Ｕ麄€計算流域范圍確定為寬8 m，高3 m，長度20 m的空間區(qū)域，通過設(shè)置對稱邊界條件模擬整個大氣環(huán)境并有效減小數(shù)值風(fēng)洞的網(wǎng)格數(shù)量。屋面設(shè)置為動網(wǎng)格邊界，墻面和地面設(shè)置為固定網(wǎng)格邊界。模型位于距離入口約6 m，距離數(shù)值風(fēng)洞出口約14 m。計算網(wǎng)格如圖7所示。

圖7 CFD計算網(wǎng)格Fig.7 CFD calculation grid

非穩(wěn)態(tài)計算選用大渦模擬。采用隨機(jī)過程模擬方法得到脈動風(fēng)速，通過將風(fēng)速作無散度化處理保證連續(xù)方程的成立，進(jìn)而作為入口風(fēng)速。本文采用文獻(xiàn)［15］中的諧波疊加法，在數(shù)值風(fēng)洞入口邊界的每一單元中心位置處進(jìn)行風(fēng)速時程模擬，程序中引入FFT算法以提高計算效率。

用上述模擬得到的時間序列還不能直接代入到實際的風(fēng)壓進(jìn)行計算，這是因為這些脈動風(fēng)速的時間序列并不一定能保證連續(xù)性方程的成立，因此，必須將邊界條件中的風(fēng)速作無散度化處理，即

由泊松方程中得到的壓力梯度修正項ΔPn代入到上述方程得到修正后的速度u(S)n+1i，為了避免流量的不平衡影響計算的收斂，在實際計算中還使用了流量修正的方法，即在計算的每個時間步上對入口處的總流量進(jìn)行修正，讓流動的進(jìn)出口流量保持一致。修正速度的大小為:

式中:ui是入口處每個單元上的脈動速度，Ai是該單元所代表的面積，A是風(fēng)速入口處的總面積，S'是該時間步出口處的流量。入口處每一單元上的脈動速度經(jīng)過無散度化后還需疊加上修正速度 Δv，這樣最后得到的脈動風(fēng)速時程才可以代入到風(fēng)荷載數(shù)值模擬計算的邊界條件上去。

在入口邊界的每個網(wǎng)格中心點處輸入一條模擬得到的風(fēng)速時程。同時將上述風(fēng)速數(shù)據(jù)經(jīng)過湍流度和流量修正，得到真正實用的非穩(wěn)態(tài)計算入口邊界。圖8所示為4 m高處一網(wǎng)格中心點位置經(jīng)過湍流度和流量修正前后的風(fēng)速時程和脈動風(fēng)速譜的比較。從圖中可以看出，經(jīng)過修正，尤其是經(jīng)過湍流度修正高處的脈動風(fēng)能量出現(xiàn)明顯下降。

圖8 入口處某點的風(fēng)速模擬與修正Fig.8 Wind velocity simulation and correction on a point at the entrance

4.2.2 有限元模型

有限元模型較為簡單，無需建立墻體等其他結(jié)構(gòu)，而僅需建立屋面模型。有限元模型前5階模態(tài)與風(fēng)洞試驗?zāi)Ｐ突疚呛希绫?所示。

表1 模型前5階自振頻率(單位:Hz)Tab.1 Top 5-order natural frequencies of the model

圖9 屋蓋的前兩階模態(tài)Fig.9 First two modes of building roof

4.3 計算結(jié)果分析

如前文所述，經(jīng)驗?zāi)B(tài)法可以提取時程信號前若干階本征模函數(shù)C1(t)～Cn(t)，其中第一個本征模函數(shù)C1(t)是從時程中分解出來的振幅最大、頻率最高的波動，依次下去的各內(nèi)在模函數(shù)，振幅逐漸變小、頻率逐漸變低。

本文提取了測點4在封閉和開孔兩種工況下下在8 m/s，9 m/s，10 m/s，11 m/s，12 m/s 五種不同風(fēng)速下的加速度時程響應(yīng)數(shù)據(jù)，聯(lián)合運用HHT和RDT法，可得到加速度時程的自由衰減曲線。對于封閉模型而言，第一、二階本征模函數(shù)還沒能把主要信號提取出來，而第三階則成功提取了響應(yīng)的典型頻率。對于開孔模型，第五階本征模函數(shù)提取了響應(yīng)的典型頻率，而從第六階開始，波動的振幅很小、頻率極低，可能是由于數(shù)據(jù)采樣頻率不夠高和波動銜接等原因造成的噪聲信號。限于篇幅，在此不列出所有自由衰減時程曲線，而由典型本征模函數(shù)運用RDT法得到的加速度時程自由衰減曲線如圖10所示。

圖10 8 m/s風(fēng)速下測點4的自由衰減曲線Fig.10 Point4 free attenuation curve at wind speed 8 m/s

對圖10所示的加速度時程自由衰減信號進(jìn)行Hilbert變換，可得到圖11所示的加速度時程振幅的對數(shù)值lna的原始曲線。

圖11 8 m/s風(fēng)速下測點4的lna曲線Fig.11 Point4 lan curve at wind speed 8 m/s

在對lna進(jìn)行線性擬合后，可計算出兩種計算工況在不同風(fēng)速下屋蓋結(jié)構(gòu)的自然頻率f和阻尼比ξ，結(jié)果見表2。

表2 各風(fēng)速下屋蓋自然頻率和阻尼比Tab.2 Natural frequency and damping ratio of roof under different wind speed

從以上圖表可以看出:① 數(shù)值計算結(jié)果的自然頻率要高于風(fēng)洞試驗，原因在于數(shù)值風(fēng)洞模型是理想密封的(即使對于迎風(fēng)面開孔工況，四周墻體是完全密閉的)，而試驗?zāi)Ｐ涂偞嬖谥欢ǖ臍怏w泄漏，結(jié)果導(dǎo)致結(jié)構(gòu)氣承剛度的降低。② 數(shù)值計算和風(fēng)洞試驗結(jié)果雖然在時程曲線上的吻合度不是非常準(zhǔn)確，但lna曲線的斜率即阻尼比的結(jié)果還是吻合得很好，這說明采用本文提出的流固耦合計算方法可以較準(zhǔn)確地計算出結(jié)構(gòu)風(fēng)振阻尼。③ 開孔狀態(tài)各個風(fēng)速下的阻尼要比封閉狀態(tài)時的大，這是由于開孔工況屋面振幅更大，從而使得風(fēng)振的氣動阻尼也更大。因此對于發(fā)生大變形的結(jié)構(gòu)振動而言，氣動阻尼的影響不可忽略。④ 由于氣動阻尼比ξa與結(jié)構(gòu)阻尼比ξs之和為總的阻尼比ξ，所以，結(jié)構(gòu)阻尼比ξs和總的阻尼比ξ已知，就可以得到氣動阻尼比 ξa，ξa=ξ-ξs。

5 結(jié)論

(1)根據(jù)三角形重心坐標(biāo)系原理提出了基于松耦合計算的流固耦合網(wǎng)格協(xié)調(diào)插值新方法，提高了流固耦合的計算效率。

(2)為解決由于柔性結(jié)構(gòu)加速度過大而導(dǎo)致的求解不穩(wěn)定問題，通過在1個時間步上進(jìn)行流體及結(jié)構(gòu)間反復(fù)交替求解來得到穩(wěn)定的時程解。

(3)HHT變換現(xiàn)已被證明是一種很好的方法而運用于結(jié)構(gòu)氣動阻尼的確定。本文通過HHT變換，很好地將振動的干擾模態(tài)分離出來，并結(jié)合RDT法，得到了理想的自由衰減曲線。

(4)計算結(jié)果顯示四周封閉建筑的阻尼比遠(yuǎn)小于迎風(fēng)面開孔建筑的阻尼比，且四周封閉狀態(tài)下氣動阻尼隨風(fēng)速變化小，其對總阻尼的影響較小;而迎風(fēng)面開孔狀態(tài)下氣動阻尼較大且隨風(fēng)速的增大而增大。因此，對于柔性且阻尼較小的結(jié)構(gòu)，特別對于開孔結(jié)構(gòu)，氣動阻尼的影響不可忽略。

［1］盧 旦，李承銘，王國儉.上海金茂大廈風(fēng)振響應(yīng)的CFD非穩(wěn)態(tài)數(shù)值模擬分析［J］.土木工程學(xué)報，2008，41(8):31-35.

LU Dan，LI Cheng-ming，WANG Guo-jian，CFD unsteady state numerical simulation of wind-induced vibration of the Shanghai Jinmao building［J］.China Civil Engineering Journal，2008，41(8):31-35.

［2］日本建筑學(xué)會.建筑風(fēng)荷載流體計算指南［M］.北京:中國建筑工業(yè)出版社，2010，12.

［3］ Gluck M，Breuer M，Durst F，et al.Computation of fluidstructure interaction on lightweight structures［J］.Wind Engineering Industry Aerodynamics，2011，89(14-15):1351-1368.

［4］Hubner B，Walhorn E，Dinkler D.Simultaneous solution to the interaction of wind flow and light weight membrane structures［C］//Bathe K J.Proceedings of the International Conference On Lightweight Structures in Civil Engineering.Warsaw，2002:519-2523.

［5］王 彬.流固耦合作用的弱耦合算法及風(fēng)與薄膜結(jié)構(gòu)的耦合分析［D］.北京:北京交通大學(xué)，2008.

［6］Batina J T.Unsteady euler airfoil solution using unstructured dynamic meshes ［J］. AIAA Journal，1990，28(8):1381-1388.

［7］ Farhat C，Degand C，Koobus B.Torsional springs for two-dimensionaldynamic unstructured fluid meshes ［J］.Computer Methods in Applied Mechanics and Engineering，1998，163(1-4):231-245.

［8］劉占生，張云峰.拉伸彈簧與扭轉(zhuǎn)彈簧模擬動態(tài)網(wǎng)格方法分析與改進(jìn)［J］.哈爾濱工業(yè)大學(xué)學(xué)報，2005，37(8):1098-1102.

LIU Zhan-shen， ZHANG Yun-feng. Analysis and enhancement on linear spring and torsional spring based dynamic mesh method［J］.Journal of Harbin Institute of technology，2005，37(8):1098-1102.

［9］Baum J，Luo H，Lohner R.A new ALE adaptive unstructured methodology for the simulation of moving bodies［J］.AIAA-94-0414，1994，32(3):414-422.

［10］ Huang N E，Shen Z，Long S R.The empirical mode decomposition and the Hilbert spectrum for nonlinear and nonstation time series analysis［J］.Proc.R.Soc.Load.A，1998，454:899-995.

［11］ Huang N E，Shen Z，Long S R.A new view of nonlinear water waves:the Hilbert spectrum ［J］.Ann.Rev.Fluid Mech.，1999，31:3417-3457.

［12］ Ibrahim SR. Random decrementtechniqueformodal identification of structures［J］.J Spacecraft，1977，14(11):696-698.

［13］劉瑞巖，張健保.隨機(jī)減量模態(tài)識別的試驗研究［J］.振動與沖擊，1993，45(1):14-9.

LIU Rui-yan，ZHANG Jian-bao.Experiment research of modal identification by random decrement［J］.Journal of Vibration and Shock，1993，45(1):14-19.

［14］楊 毅.屋蓋結(jié)構(gòu)的風(fēng)振響應(yīng)與氣動阻尼［D］.杭州:浙江大學(xué)，2004.

［15］王之宏.風(fēng)荷載的模擬研究［J］.建筑結(jié)構(gòu)學(xué)報，1994，15(1):44-52.

WANG Zhi-hong.Simulation of wind loading［J］.Journal of building structures，1994，15(1)，44-52.