基于大渦模擬的臺北101大樓風致響應分析

盧春玲， 李中洋， 李秋勝

(1.桂林理工大學 土木與建筑工程學院，桂林 541004；2.桂林理工大學 廣西有色金屬隱伏礦床勘查及材料開發協同創新中心，桂林 541004； 3.廣西巖土力學與工程重點實驗室，桂林 541004； 4.香港城市大學 建筑系 香港 100013)

1 工程概況

超高層建筑具有質量輕、強度高、高寬比大且自然振動頻率以及阻尼較低等特性。故對風力所造成的擾動較敏感，在風力作用下可能產生較大的風致響應。基于安全性與使用者的舒適性考慮，風力設計往往成為決定超高層建筑結構設計的重要因素。

臺北101大樓地面以上101樓，地下5層。地面以上高度508 m，坐落在臺北繁華地段信義計劃區內。

大樓的立面圖如圖1所示，大樓造型宛若勁竹節節高升。平面為正方形，邊長從46～64 m不等。大樓采用雙層結構，外層由8根巨大的鋼柱組成，在大樓的四個外側分設兩根。巨型柱的最大截面長3 m、寬2.4 m，自地下5樓貫通至地上90樓，柱內灌入高密度混凝土，外以鋼板包覆。地震力和風力作用時，這些巨型柱提高大樓的橫向剛度。內層結構為大樓提供可用空間。大樓每隔8到10層設置鋼桁架組成的轉換層。圖2給出了桁架轉換層的位置。

圖1 臺北101大樓立面圖

圖2 臺北101大樓結構體系立面圖

為了降低大樓受高空強風及臺風作用造成的搖晃，大樓內設置了單擺式TMD(Tuned Mass Damper, TMD)系統，即調頻質量阻尼器。它由數塊鋼板堆疊焊接而成。鋼球直徑為5.5 m，總重約為660 t，鋼球底下設置8根液流阻尼器，以提供TMD系統額外的消能機制，以降低鋼球晃動行程。質塊阻尼器利用鋼索懸吊，并由其調整鐘擺的擺長(TMD懸吊于92樓至88樓)，使其自振頻率與結構的基本頻率一致而達到吸收結構振動能量的功能。阻尼器的設置詳見圖3。

臺北101大樓屬風敏感建筑。同時臺北市又位于世界上臺風多發的區域，臺風可能使得101大樓承受遠遠大于普通高層建筑承受的荷載。這些特點使得該建筑需要對臺風作用下大樓的結構性能進行深入細致的研究。

(a)總體圖(b)詳細構造圖

圖3 臺北101大樓阻尼器布置圖

Fig.3 Damper layout of Taipei 101 Tower

大渦模擬(LES)是現今計算風工程的研究熱點之一。隨著數值模擬技術和計算機性能的提升，研究人員嘗試將大渦模擬運用到工程實踐中。例如，顧明[1-3]、李秋勝[4-5]、周志勇[6]等采用LES對結構風效應進行了分析[7-9]。利用大渦模擬進行實際超高層建筑中進行結構風荷載及響應研究還較少[10]。本文應用大渦模擬并結合一種新的可滿足大氣邊界層中風場特性的湍流脈動速度生成方法——離散再合成的隨機湍流生成法(DSRFG)[11]模擬非穩態邊界層湍流風場。在大渦模擬的亞格子模型方面，采用一種新亞格子模型[12]，基于Linux系統下Fluent軟件，模擬得到臺北101大樓周圍的風流場及作用于其上的風荷載時程數據。建立大樓的3維有限元計算模型，計算出風荷載作用下大樓的風致響應以及等效靜力風荷載。并將計算結果與現場實測以及風洞試驗的相應數據進行對比，以對該數值模擬方法的實用性和準確性進行檢驗。

2 大樓強振監控系統概況

安裝在臺北101大樓上的強振監測系統主要由傳感器和數據采集處理系統兩部分組成。傳感器系統包含30個加速度傳感器，分別安裝在大樓的-5th (地下室)，1 st，36 th，60 th，89 th和101 st層，如圖4所示。為測得大樓的扭轉加速度，部分加速度傳感器設置在所在層的中心為對稱點的兩個不同位置，如圖5所示。通過分析和處理傳感器測得的數據，可獲得大樓的加速度響應以及動力特性，如自振頻率、模態以及阻尼比等數據。

本文用于對比分析的現場實測數據是由采樣頻率為200 Hz的數據采集系統從2005年8月到2008年5月期間，三次臺風作用和一次地震作用下采集得到。

圖4 加速度傳感器立面布置圖

圖5 加速度傳感器平面布置圖

3 大樓大渦模擬

3.1 模型建立與網格劃分

臺北101大樓的計算模型為全尺寸，數值模擬的計算模型以及計算域如圖6所示，X方向長為23 Db(Db為建筑的寬度)，Y方向寬為16 Db，Z方向高度為兩倍的建筑高度。計算域中建筑物的阻塞比小于3%。

網格的劃分情況如圖7所示，對建筑物周圍網格進行局部加密，為了提高網格質量，這部分網格采用四面體網格。在加密區域以外采用了六面體網格。總的網格數量為320萬左右。建筑模型表面采用5層棱柱體網格，第一層壁面網格的高度選取為0.01 m，最大網格歪斜度(skewness of equisize angle)為0.89，網格尺寸增長率取為1.05。近墻壁y+值約為35～85，適用于壁面函數。

圖6 臺北101大樓計算模型及計算域

圖7 臺北101大樓計算模型及計算域

3.2 邊界條件的設定

本文大渦模擬采用作者提出的一種新的湍流入口生成方法-DSRF方法[11]。此方法在綜合了以往方法優勢的基礎上，具備的優點包括：①嚴格保證入口湍流滿足連續性條件div(u)=0；②基于嚴格的理論推導，具有通用性。生成的脈動速度滿足指定的譜密度函數；③入口湍流的空間相關性可通過相關性尺度因子調整；④每個坐標點的入口湍流生成過程相互獨立，適用于并行計算。⑤能處理輸入湍流功率譜及湍流積分尺度的各向異性。速度入口處平均風速沿高度的變化服從指數率，見式(1)。

(1)

式中：Z10和V10為10 m高度和10 m高度處的風速。基于臺北松山機場氣象站測量的近地風速觀測數據采用風氣候統計模型得到50年回歸期的十分鐘平均風速V10為每秒43.27 m/s，α=0.15，Z為高度變量。

本文基于Shiau[13]在臺灣基隆港的觀測塔測(距離臺北北部約20 km)得的26 m高度處的強風(臺風)數據的結果：縱向湍流積分強度在0.18～0.23之間，縱向湍流積分尺度在40～200 m之間，對應的風速為27～45 m/s。以及李秋勝等[14-15]在強風作用時，現場實測得到香港中環大廈(374 m)以及深圳定王大廈(384 m)樓頂處的風場數據結果：湍流強度在0.1～0.4之間，湍流積分尺度在171～600 m之間，風速在5～25 m/s。基于以上的現場實測數據，確定了臺北101大樓入口處的湍流強度剖面及湍流積分尺度的分布，如圖8所示。采用DSRFG方法生成的入口處的瞬時風速分布圖如圖9所示。通過DSRFG方法生成的入口處脈動風速時程的功率譜與目標譜卡曼譜的比較如圖10所示。由圖10可知，大渦模擬所生成的脈動風速時程功率譜與目標譜一致，故其功率譜特性符合大氣邊界層湍流的要求。

圖8 入口湍流強度及積分尺度剖面

圖9 DSRFG方法產生的入口瞬時風速場

Fig.9 Instantaneous velocity field at the inflow boundary generated by the DSRFG method

圖10 DSRFG方法生成的脈動風速時程功率譜與目標譜比較

Fig.10 Velocity sequence spectrum generated by DSRFG and comparison with the target spectrum

流場的出口認為湍流完全發展，流場變量散度為零。計算域頂部與兩側為自由滑移壁面，壁面剪應力為零。建筑表面和地面采用無滑移的壁面條件限定流體和固體區域。

3.3 大渦模擬亞格子模型及計算求解

本文采用的大渦模擬新的亞格子模型的詳細介紹參見文獻[12]。新亞格子模型的主要特點是：適合工程應用、一方程模型、不需采用試驗濾波, 適合低階格式及無結構網格應用，具有動態特征，且計算量少。

在數值求解中，認為空氣的流動是不可壓縮的。數值計算采用壓力隱式分割算法(PressureImplicit with Splitting of Operators, PISO)進行迭代求解。時間離散采用二階隱式格式，空間離散采用二階中心格式。為了使得計算更快的收斂，在大渦模擬計算之前先進行了RANS模型的計算，將RANS模型計算的結果通過瞬態化處理作為大渦模擬計算的初始流場。考慮到大渦模擬計算的準確性及計算時間成本，大渦模擬非定常計算的時間步長取為0.05 s。共進行了16 000步的大渦模擬非定常計算。流場數據的統計計算采用了最后11 000步的計算結果。

本文的計算在32CPUs并聯成的高性能計算機集群上進行。LES的計算開銷如下：內存使用12 GB，計算總耗時達733 h，每個時間步長耗時約60 s。速度和壓力收斂殘差標準選取為10-5。

4 大渦模擬計算結果與分析

4.1 流場

由于風洞試驗受到試驗條件和手段的限制，很難測得建筑物周圍的整個流場的情況，CFD方法則不受此限制，可以給出非常完整的結果。如風速分布和風速矢量圖等，為設計提供參考。

大渦模擬得到的臺北101大樓周圍的平均、瞬時以及根方差風流場如圖11(a)～(c)所示。從圖中可以看出，典型鈍體擾流的流場特點，如在建筑物迎風前方以及兩側的流動分離、渦旋脫落、尾流等都得到真實的再現。從圖11(b)可以看出，在建筑物前方有大量的隨機渦旋結構。不同頻率的風速脈動與分離區域中的渦旋結構相互作用，產生密集的渦旋脫落，這也是高層建筑橫風向脈動風力產生的主要原因。

根方差風速代表了風速脈動大小，從圖11(c)中可以看出，建筑物前方的根方差風速較小，其主要由來流湍流產生。在來流遇到建筑物的阻擋后，由于建筑物對流場的擾動，產生特征湍流，較大的根方差風速主要出現在流動分離區域以及渦旋脫落的尾流區域。

4.2 風荷載

將臺北101大樓分51層(基本兩個樓層為一層)，通過積分得到這51層處大樓兩個主軸方向上的風荷載Fx，Fy以及繞Z軸的扭矩的時程數據。通過對51層處的風荷載數據進行積分還可得到大樓總的風荷載。臺北101大樓總的順風向、橫風向以及扭矩風荷載時程曲線見圖12。從圖12可看出風荷載的脈動是無序的，這說明數值風洞模擬預測出了作用在大樓上的隨機風力。

(a)平均風速流場

(b)瞬時風速流場

(c)根方差風速流場

圖11 臺北101大樓周圍風速流場

Fig.11 Instantaneous velocity contours around Taipei 101 Tower

(a)順風向風荷載

(b)橫風向風荷載

(c)扭轉向風荷載

圖12 臺北101大樓總的風荷載時程曲線

Fig.12 Time-history of the wind forces

臺北101大樓總的風荷載的功率譜以及2008年3月12日“汶川”地震作用下大樓實測得到的基底加速度功率譜如圖13及14所示。對比圖13和圖14發現：大樓總的風荷載的功率譜能量與遠場地震波導致的樓底的加速度的能量都主要集中在0～2 Hz頻段。造成這一特殊現象的原因可能是地震波的高頻部分在遠距離(約1 890 km)的地層傳輸過程中被過濾掉了。這一現象也說明遠場地震波對大樓的激勵在一定程度上跟風荷載對大樓的激勵相類似。

(a)順風向風荷載

(b)橫風向風荷載

(c)扭轉向風荷載

圖13 臺北101大樓總的風荷載功率譜

Fig.13 Power spectral densities of the overall wind forces

圖14 “汶川”地震作用下臺北101大樓基底加速度功率譜

Fig.14 Power spectral densities of the accelerations at the basement (B5) during Wenchuan earthquake

大樓在50層、40層、30層以及20層處積分得到的順風向和橫風向以及扭轉的風荷載的功率譜分別如圖15～17所示。

從圖15中看出：順風向風荷載的功率譜與典型的順風向風荷載功率譜類似，能量集中在較寬的頻帶范圍。不同層上的順風向風荷載功率譜的峰值隨著樓層的增高而下降，也就是低層的順風向風荷載功率譜的最大能量高于高層的。這可能是由于氣流在建筑的迎風面駐點以下下沉至地面形成回流，這一區域的湍流較大所造成的。除了功率譜峰值頻段，不同層上的風荷載功率譜大小沿樓高沒有明顯的變化。

從圖16表示的各層上橫風向風荷載功率譜曲線來看，每條功率譜曲線都存在一個明顯的尖峰，這可能是渦旋脫離產生的。這也說明橫風向風荷載主要是渦旋脫落引起的。不同層上橫風向功率譜的峰值沿建筑高度沒有明顯的變化。然而，在較高頻段，各條功率譜曲線隨著層高增加而抬高。這一現象的合理解釋為：在較高的樓層，3維流動效應以及較高的風速在大樓側墻上產生了較之低處較高的風壓脈動。這也是橫風向風荷載功率譜中高頻部分能量的主要產生原因。

大樓各層上的扭矩功率譜如圖17所示，在各條功率譜曲線上也都有一個明顯的峰值，這一峰值所對應的頻率也跟橫風向荷載功率譜的峰值頻率比較接近，這說明由來流風湍流和渦旋脫落導致的不平衡的脈動風壓是產生扭矩的重要原因。

圖15 臺北101大樓總的順風向風荷載功率譜

圖16 臺北101大樓總的橫風向風荷載功率譜

圖17 臺北101大樓總的扭轉風荷載功率譜

5 風致響應分析

5.1 大樓三維有限元模型

臺北101大樓的三維有限元模型如圖18所示。在建立有限元模型時采用了4類單元：3維梁單元模擬柱和梁。3維桿單元模擬鋼支撐。樓板用板單元進行模擬。活載以及其他非結構部分采用3維質量單元模擬。大樓在基礎處的約束為3向固定約束。大樓的三維有限元模型共計20 532個梁單元，24 048個板單元, 以及3 496個桿單元。

5.2 模態分析

將有限元模型模態分析計算的大樓X和Y方向前三階以及第一階扭轉自振頻率與臺風“馬莎”、“泰利”、“柯羅莎”及“汶川”地震作用時臺北101大樓現場實測數據得出的自振頻率的平均值進行對比，結果列于表1中。從表1中可看出：有限元計算的自振頻率與實測結果相差在5.9%～14%之間。這表明可利用有限元模型對大樓的風致響應進行分析。

基于有限元模型根據模態分析得到的臺北101大樓前3階X和Y方向自振頻率對應的振型以及第一階扭轉自振頻率對應的振型(相對于重心軸)如圖19所示。

圖20以及圖21給出了臺北101大樓由有限元模型計算出的以及通過3次臺風以及“汶川”地震時由實測加速度數據得到的前三階X、Y方向的振型(取大樓第101樓層的振幅為1)。從圖中可看出有限元分析得到的振型與現場實測的結果基本一致。

圖18 有限元模型

5.3 位移和加速度響應

通過ANSYS有限元瞬態動力分析得到的順風向(X軸方向)以及橫風向(Y軸方向)臺北101大樓最高的居住層(第89層)中心處的位移以及加速度響應的時程序列如圖22～23所示。

一階 二階 三階(a)X向一階 二階 三階(b)Y向一階(c)扭轉

圖19 臺北101大樓計算模態

(a)3次臺風和“汶川”地震中實測X向振型

(b)有限元計算的X向振型

(a)3次臺風和“汶川”地震中實測X向振型

(b)有限元計算的X向振型

風荷載作用下，臺北101大樓每層樓中心處X及Y軸方向上的位移以及加速度響應如圖24～25所示。從圖中可以看出：①大樓風致響應曲線較光滑，沒有明顯的拐點，其曲線形狀類似于均布力作用下懸臂梁的變形曲線。②橫風向的最大以及根方差位移以及加速度響應值均大于順風向的結果。這也說明臺北101大樓的橫風向響應起主導作用。因此，在大樓的抗風設計中應該考慮橫風向的響應以及風荷載。

(a)順風向

(b)橫風向

(a)順風向

(b)橫風向

表2給出了風荷載作用下，計算得到的臺北101大樓最高居住樓層(第89層，離地面高度382.2 m)中心處最大加速度響應的結果，同時，風洞試驗的結果也列于表2。從表2可以看出：數值計算的結果與風洞試驗的結果比較一致，計算結果稍大于試驗結果。兩者之間的相對誤差相差小于13%。

(a)順風向

(b)橫風向

圖25 大樓立面上加速度曲線

計算值/(m·s-2)風洞試驗/(m·s-2)差別/%阻尼比風速重現期 (年)1.5%5%1.130.651.0030.5812.61210

為了研究超高層建筑的風致響應的特性，對大樓不同樓層的位移以及加速度響應進行了功率譜分析。不同樓層處X及Y軸方向上的位移以及加速度響應功率譜如圖26～27所示。從圖中可看出在0.17 Hz處，各層樓上的位移以及加速度功率譜出現一明顯的峰值，這一頻率恰好對應大樓的X及Y軸方向上的一階自振頻率。不同層上的位移及加速度功率譜的峰值隨著樓層的增高而增加。除了0.17 Hz出現峰值，在0.39 Hz以及0.41 Hz處，順風向及橫風向加速度功率譜還出現了第二個較小的峰值，這兩頻率也與大樓在這兩個方向上的二階自振頻率一致。由此可見：大樓的X及Y軸方向上的位移、加速度響應以一階頻率為主，而二階頻率對大樓加速度響應也有一定影響。

(a)X方向

(b)Y方向

(a)X方向

(b)Y方向

計算與現場實測得到的大樓最高居住層第89層處，X及Y軸方向上的加速度功率譜對比如圖28所示。從圖28中可看出，計算的加速度功率譜與實測譜吻合較好。但在高頻部分，計算的加速度功率譜小于實測譜，造成這一差別的原因可能是數值模擬中沒有考慮周邊建筑，實際風場的湍流強度與數值模擬存在一定差別。

(a)X方向

(b)Y方向

Fig.28 Power spectra of the measured and calculated accelerations at the top occupied floor

6 結 論

(1)入口湍流產生的DSRFG方法可以產生滿足大氣邊界層風場特性的卡曼功率譜、LES入口的順風向湍流強度剖面以及順風向湍流積分尺度豎向分布的湍流風場。

(2)基于大渦模擬得到的風荷載時程數據作用于臺北101大樓的3維有限元模型經過瞬態動力計算得到的大樓風致響應結果與現場實測以及風洞試驗結果吻合較好。說明本文的數值模擬方法可以較為準確地模擬高雷諾數流場。

(3)總的來說，運用新的大渦模擬亞格子模型結合合適的湍流入口生成方法能夠對高層建筑的風荷載及風致響應進行有效的預測。數值模擬計算結果通過與現場實測與風洞試驗結果的對比證明其有較好準確性。數值模擬可為高層建筑結構抗風設計提供有效參考。