基于大渦模擬超大直筒-錐段型鋼結構冷卻塔極值風壓研究

王 浩, 柯世堂

(南京航空航天大學 航空宇航學院， 南京 210016)

作為一種新穎的冷卻塔結構形式，大型鋼結構冷卻塔屬于典型的高聳空間網架結構，較混凝土冷卻塔自重更輕、柔性更大、阻尼更小，風荷載已成為此類冷卻塔結構設計中的控制荷載。同時，直筒-錐段型鋼結構冷卻塔在用鋼量、加工制造和施工安裝等方面較雙曲線型鋼結構冷卻塔具有優勢，在國際上已有十余例成功應用于電廠間接空冷建設，而在國內此類鋼結構冷卻塔尚處于預研階段[1]。

對于直筒-錐段型鋼結構冷卻塔這類典型風敏感結構，研究人員和工程師們最關注的往往是風荷載的極值，尤其是風壓極值作用下冷卻塔的局部和整體穩定性能驗算[2-4]。而現有冷卻塔設計規范[5]和建筑結構荷載規范[6]中僅給出了此類圓柱形斷面沿環向的平均風壓分布曲線；柯世堂等[7]基于改進的峰值因子法探討了大型雙曲冷卻塔表面極值風壓的取值問題；李鵬飛等[8]在考慮相關性的基礎上給出了雙曲線冷卻塔外表面極值風壓建議取值；羅堯治等[9]基于計算流體動力學(Computational Fluid Dynamics, CFD)方法進行了屋蓋結構的風壓極值和相關性分析。由于直筒-錐段型鋼結構冷卻塔特殊的結構形式和氣動外形，導致難以尋找相應的設計參考標準，現有規范和文獻均沒有明確給出此類冷卻塔的風壓極值。

鑒此，以國內擬建的首座超大直筒-錐段型鋼結構冷卻塔為例，運用大渦模擬(Large Eddy Simulation, LES)法模擬出該塔外表面三維氣動力時程；并通過與國內外大型冷卻塔實測及風洞試驗結果對比驗證了數值模擬結果的有效性。在此基礎上，基于三種常用極值估算方法(峰值因子法[10]、改進峰值因子法[11]和Sadek-Simiu法[12])對比分析了此類鋼結構冷卻塔表面峰值因子和風壓極值的分布規律，并基于最小二乘法分別給出了直筒和錐段部分的風壓極值一維/二維擬合公式及取值建議。相關研究結論可為此類直筒-錐段型鋼結構冷卻塔的抗風設計提供參考依據。

1 工程概況與大渦模擬

國內某擬建的直筒-錐段型鋼結構冷卻塔塔高189 m，由主筒、內部加強桁架、外部附屬桁架和外表面的圍護層組成。進風口高32.5 m，直徑123.2 m，進風口高度以上結構外覆擋風鋼板，表面光滑平整。直筒-錐段型鋼結構冷卻塔主要結構尺寸與構件信息，如表1所示。

表1 直筒-錐段型鋼結構冷卻塔主要結構尺寸與構件信息表

該冷卻塔幾何模型按實際尺寸建立，以30%透風率[13]考慮正常運行狀態下的百葉窗開啟效應。塔筒外表面沿環向逆時針每間隔12°布置一個測點，沿子午向不同高度布置16層測點，共計16×30=480個測點。為保證尾流充分發展，數值風洞尺寸定為24D×15D×4H(順風向X×橫風向Y×高度方向Z，D為底部直徑，H為塔高)，冷卻塔模型位置距離出口17D，阻塞率為1.21%。將整個計算域劃分為外圍區域和局部加密區域以兼顧計算效率和精度，外圍區域形狀規整，采用高質量的結構化網格進行劃分，局部加密區域內含冷卻塔模型，采用非結構化網格進行劃分。加密區最小網格尺寸為0.5 m，總網格數量約640萬個，網格質量>0.4。冷卻塔測點布置及網格劃分，如圖1所示。

定義進口邊界條件為速度入口，按照B類地貌設置大氣邊界層指數風速剖面和湍流度剖面，地面粗糙指數為0.15，10 m高度處的基本風速為26.83 m/s，通過用戶自定義函數(User Defined Function)實現上述入流邊界條件與FLUENT的連接，如圖2所示。出口采用壓力出口邊界條件；頂部和側面采用對稱邊界條件；地面及冷卻塔表面采用無滑移壁面，如圖3所示。基于LES法對冷卻塔周圍復雜的流場進行模擬。亞格子尺度選用Smagorinsky-Lilly模型，計算殘差設置為1×10-6，時間步長定為0.05 s。

2 大渦模擬結果與有效性驗證

2.1 數值模擬方法有效性驗證

由于現階段尚無直筒-錐段型冷卻塔的風荷載研究見于文獻，為驗證本文數值模擬方法的有效性，首先進行了相同塔高的雙曲線型冷卻塔風荷載CFD模擬。圖4給出了該塔塔筒喉部外表面體型系數的平均值和根方差沿環向分布曲線，并分別在圖中與水工規范[10]和大型雙曲冷卻塔實測及風洞試驗結果[14-15]進行對比。由圖4可知，數值模擬得到體型系數平均值與規范曲線吻合良好，且根方差分布曲線與Sageau進行的實測得到的曲線較為接近，考慮到風壓信號的脈動程度與實測塔所處的地形、來流湍流和周邊干擾密切相關，因此可認為LES法具有一定的有效性。

(a) 測點布置

(b) 整體網格

(c) 局部網格

圖2 風速剖面及湍流強度分布

圖3 計算區域與邊界條件示意圖

2.2 直筒-錐段型冷卻塔

前文雙曲線冷卻塔風荷載大渦模擬結果準確可信的基礎上，采用相同的參數設置進行了直筒-錐段型鋼結構冷卻塔風荷載數值模擬。直筒-錐段型冷卻塔典型截面表面繞流及渦量分布圖，如圖5所示。由圖5可知，氣流在冷卻塔迎風面發生流動分離且局部出現加速效應，氣流流動分離后在冷卻塔頂部和側面產生漩渦脫落，持續發展后在冷卻塔背風面形成尾流渦旋及回流。

(a) 平均值驗證

(b) 根方差驗證

分別對直筒和錐段外表面各層體型系數平均值和根方差進行平均，并與水工規范、建筑結構荷載規范和國內外大型雙曲冷卻塔實測及風洞試驗結果進行對比，如圖6所示。對比分析可知：① 錐段部分體型系數分布曲線的負壓極值點和分離點對應角度與水工規范中雙曲線冷卻塔風壓分布曲線較為一致，直筒段體型系數分布曲線與荷載規范中圓截面構筑物分布曲線基本吻合；② 脈動風壓沿環向分布規律與國外雙曲冷卻塔的實測結果較為接近，但在側風面和背風面的數值要稍大于實測結果，考慮到直筒-錐段型鋼結構冷卻塔氣動外形與雙曲冷卻塔差異明顯，本文基于大渦模擬得到的脈動風壓具有一定的可信度，可用于后續的極值風壓研究。

(a) 直筒頂部流線圖

(b) 側立面流線圖

(c) 直筒頂部渦量圖

(d) 側立面渦量圖

需要說明的是，本文后續研究給出的直筒和錐段峰值因子和極值風壓等數據均是對直筒部分和錐段部分所有測點層進行平均的結果。

(a) 平均值驗證

(b) 根方差驗證

3 極值估計方法

工程設計中脈動風荷載和風振響應極值都應考慮一定的保證系數，計算公式如下

*g*σr

(1)

3.1 峰值因子法

早期Davenport等以高斯分布假定為基礎，利用零值穿越理論給出一峰值因子，以該峰值因子乘以風壓的根方差后再加上平均風壓得到極值風壓，此即為峰值因子法。峰值因子按下式計算

(2)

式中:r=0.577 2，為歐拉常數；v表示單位時間內數據穿越平均值的次數；T為計算的時間長度。

這類方法為不少國家的荷載規范所引用，但對于時間T的取值方式略有不同，我國規范的風壓時距均為10 min。按照此方法得到的峰值因子偏小，約為3.0～3.5。

3.2 改進峰值因子法

Winterstein等[16]在Davenport的工作基礎上，將高斯隨機變量表示為考慮高階統計量(偏度和峰度)的非高斯隨機變量的Hermite多項式，從而將Davenport僅僅適用于高斯過程的基于零值穿越理論的峰值因子法擴展到適用于估算非高斯過程的風壓極值，此即為峰值因子法。給出非高斯隨機變量x(t)和標準高斯過程u(t)的關系式

x=α{u+h3(u2-1)+h4(u3-3u)}

(3)

通過Hermite多項式，可得非高斯過程的峰值因子

(4)

3.3 Sadek-Simiu法

隨著統計學的發展以及在工程上的應用，Sadek和Simiu在零值穿越理論的基礎上提出了非高斯過程的極值計算方法(Sadek-Simiu法)。時距T內，時程y(t)的極值ypk,T的概率分布函數為

(5)

(6)

式中：v0,y為高斯過程y(t)的零值穿越率，與峰值因子法的含義相同。

因為在映射過程中并不能得到實際的高斯時程y(t)，因此上式在計算時仍然由x(t)的譜Sx(n)代入計算。峰值累積分布函數Fyoi,r(ypi,T)可由式(5)得到，將高斯過程映射到非高斯過程上即可獲得非高斯分布時間歷程x(t)的峰值累積分布。

4 極值風壓研究

4.1 非高斯特性

通過繪制頻率分布直方圖獲得測點風壓概率密度分布，并與標準高斯分布曲線進行比較發現環向斷面上風壓既有高斯分布又包含非高斯分布。限于篇幅，本文僅給出環向典型非高斯測點的風壓時程曲線和概率密度分布曲線，如圖7和圖8所示。由圖可知，冷卻塔側風面位于來流風分離區，風壓時程曲線具有明顯的不對稱性，伴隨有大量間歇性脈沖信號并伴隨有很大的風吸力。從更多測點體型系數的概率分布看，結構表面有大部分區域的風壓信號概率分布屬于非高斯分布。

(a) 負壓極值區

(b) 分離區

(a) 負壓極值區

(b) 分離區

4.2 峰值因子取值探討

由三種極值風壓算法計算得出的直筒段和錐段沿環向峰值因子平均值，如圖9所示。對比可以發現三種不同方法計算得出的結果差異較大，基于高斯分布假定的峰值因子法計算得出的數值普遍偏低，且沿環向角度變化微弱；錐段部分背風區的峰值因子較直筒段背風區更大，說明此類鋼結構冷卻塔在錐段區域的背風面可能會產生更為明顯的風壓極值。

(a) 錐段部分

(b) 直筒部分

不同算法計算得到的峰值因子二維分布圖，如圖10所示。對比發現峰值因子分布呈明顯的二維分布特征，改進峰值因子法和Sadek-Simiu法計算得到的峰值因子沿高度變化明顯，尤其是在錐段和直筒頂部等三維繞流明顯的區域峰值因子變化顯著。三種算法計算得到的峰值因子最大值均出現在背風區，這是由于背風區風壓時程的均值較小，且對應脈動值相對較大。

(a) 峰值因子法

(b) 改進峰值因子法

(c) Sadek-Simiu法

4.3 風壓極值對比研究

三種方法計算得出的直筒段和錐段體型系數極值曲線，并與規范曲線進行對比，如圖11所示。對比發現錐段迎風面和側風面的極值風壓接近于雙曲冷卻塔推薦曲線，而直筒段迎風面和側風面的極值風壓與大型圓截面構筑物風壓極值曲線更接近。三種方法計算得到的錐段背風區的極值風壓分別比規范數值大35.3%、58.9%和39.5%，直筒段計算結果分別比規范數值大37.4%、54.9%和45.6%。

(a) 錐段部分

(b) 直筒部分

Fig.11 The comparison diagram of extreme shape coefficient under different methods with standard value

直筒-錐段型冷卻塔體型系數極值二維分布圖，如圖12所示。由圖12可知，① 極值風壓分布關于風軸方向比較對稱，最大正值和負值分別位于冷卻塔的迎風面和側風面，最大負值集中分布于60～125 m高度區間內；② 錐段側風區風壓負值大小及分布面積明顯小于上部直筒段，背風區域極值風壓分布也與直筒段有很大區別。對比表明改進峰值因子法計算得出的風壓極值基本包絡住其它兩種方法的結果和規范值，后續極值擬合研究均基于改進峰值因子法的計算結果。

4.4 極值風壓擬合

上述研究表明直筒-錐段型鋼結構冷卻塔極值風壓沿環向和子午向分布呈現出典型的二維特征，而相關設計規范只給出了雙曲線冷卻塔和圓截面構筑物的的一維風壓極值曲線，該方法對直筒-錐段型鋼結構冷卻塔而言不盡合理。為方便工程研究與設計人員精確獲得此類冷卻塔極值風壓，本文基于非線性最小二乘法原理，以子午向高度比n和環向角度θ為目標函數，擬合給出此類冷卻塔二維風壓極值的計算公式，同時給出了以環向角度θ為目標函數的直筒段和錐段分段一維平均擬合公式，供設計參考。

錐段和和直筒段原始數據與一維擬合曲線對比示意圖，如圖13所示。由圖13可知，公式擬合結果與實際數據一致，擬合公式定義如式(7)所示。具體參數見表2和表3。

(7)

(a) 峰值因子法

(b) 改進峰值因子法

(c) Sadek-Simiu法

(a) 錐段

(b) 直筒部分

參數a1a2a3b1b2b3c1c2c3數值2．40300．0093-3．24105．96900．04200．28875．61900．04393．0860

表3 直筒部分目標擬合公式參數表

直筒-錐段型鋼結構冷卻塔二維極值風壓擬合公式定義為

Cpm(n,θ)=a00+a10θ+a01n+a20θ2+a11θn+a02n2+

a21θ2n+a12θn2+a03n3+a22θ2n2+a13θ1n3+

a04n4+a23θ2n3+a14θn4+a05n5

(8)

式中：n=z/H，其中z為測點對應高度；H為塔頂總高度；aij(i=1～2,j=1～5)為擬合系數；具體數值，見表4。

直筒-錐段型鋼結構冷卻塔極值風壓實際分布與二維擬合曲面對比圖，如圖14所示。圖14中散點值即為冷卻塔真實極值風壓，曲面對應數值為根據二維擬合公式計算得到的極值風壓。從擬合曲面的整體分布來看，其沿子午向和環向的變化規律與實際極值風壓基本一致，數值上可以包絡住真實數值，誤差率<10%。對比結果表明本文提出的二維極值風壓擬合公式可作為此類大型直筒-錐段型鋼結構冷卻塔極值風壓取值的參考依據。

表4直筒-錐段型冷卻塔二維極值風壓目標擬合公式系數表

Tab.4Twodimensionalextremewindpressurefittingformulacoefficienttableofcylindrical-conecoolingtower

擬合系數a001．245a0532．4400a140．003036a0113．61a10-0．1440a200．001526a02-64．49a110．0332a21-0．000109a03123．30a12-0．0690a220．000219a04-104．30a130．0408a23-0．000129

圖14 極值風壓原始數據與二維擬合曲面對比圖

Fig.14 Comparison between the raw data and the two dimensional fitted surface of extreme wind pressure

5 結 論

本文系統詳盡地研究了直筒-錐段型鋼結構冷卻塔表面風壓極值分布特性，主要涉及大渦模擬、峰值因子計算方法、風壓極值計算方法和非線性最小二乘法擬合公式。得到如下主要結論：

(1) 通過與國內外實測及風洞試驗結果對比驗證了LES方法和結果的有效性，風荷載時程可用于后續風壓極值研究。

(2) 直筒-錐段型鋼結構冷卻塔負壓極值區和分離區附近風荷載呈現明顯的非高斯特性，基于高斯分布假定的峰值因子法計算得出的極值風壓數值偏小，改進峰值因子法計算結果基本包絡住其他兩種方法，Sadek-Simiu法計算得到的風壓極值介于二者之間。

(3) 錐段和直筒段背風面極值風壓均大于規范數值，錐段迎風面和側風面極值風壓與水工規范極值接近，而直筒段迎風面和側風面極值風壓與荷載規范極值接近。

(4) 本文提出的二維風壓極值擬合公式可包絡直筒和錐段的實際風壓極值，最大誤差控制在<10%，可對此類鋼結構冷卻塔風荷載設計取值預測提供參考。

[1] 鄒云峰, 李壽英, 牛華偉,等. 雙曲冷卻塔等效靜力風荷載規范適應性研究[J]. 振動與沖擊, 2013, 32(11):100-105.

ZOU Yunfeng， LI Shouying, NIU Huawei, et al. The hyperbolic cooling tower equivalent static wind load code for the adaptive study[J]. Journal of Vibration and Shock, 2013, 32 (11): 100-105.

[2] 柯世堂, 葛耀君, 趙林. 大型雙曲冷卻塔表面脈動風壓隨機特性—非高斯特性研究[J]. 實驗流體力學, 2010, 24(3): 12-18.

KE Shitang, GE Yaojun, ZHAO Lin. Research on festures of fluctuating wind pressure on large hyperbolic cooling tower: features of non-Gaussian[J]. Journal of Experiments in Fluid Mechanics, 2010, 24(3): 12-18.

[3] 全涌, 顧明, 陳斌,等. 非高斯風壓的極值計算方法[J]. 力學學報, 2010, 42(3):560-566.

QUAN Yong, GU Ming, CHEN Bin, et al. Study on the extreme value estimating method of non-Gaussian wind pressure[J]. Chinese Journal of Theoretical and Applied Mechanics, 2010, 42(3): 560-566.

[4] GIOFFRE M, GUSELLA V. Damage accum ulation in glass plates[J]. Journal of Engineer Mechanics, 2002(7): 801-805.

[5] 工業循環水冷卻設計規范:GB/T 50102—2014[S]. 北京: 中國計劃出版社, 2014.

[6] 建筑結構荷載規范:GB 50009—2012[S]. 北京: 中國建筑工業出版社, 2012.

[7] 柯世堂, 葛耀君, 趙林. 大型雙曲冷卻塔表面脈動風壓隨機特性——風壓極值特性探討[J]. 實驗流體力學, 2010, 24(4): 7-12.

KE Shitang, GE Yaojun, ZHAO Lin. Research on festures of fluctuating wind pressure on large hyperbolic cooling tower: discussions on extreme wind pressure[J]. Journal of Experiments in Fluid Mechanics, 2010, 24(4): 7-12.

[8] 李鵬飛, 趙林, 葛耀君, 等. 超大型冷卻塔風荷載特性風洞試驗研究[J]. 工程力學, 2008, 25(6):60-67.

LI Pengfei, ZHAO Lin, GE Yaojun, et al. Wind tunnel investigation on wind load characteristics for super large cooling towers[J]. Engineering Mechanics, 2008, 25(6):60-67.

[9] 羅堯治, 丁慧. 波紋表面月牙形懸挑屋蓋風荷載特性[J]. 浙江大學學報(工學版), 2013, 47(6):1063-1071.

LUO Yaozhi, DING Hui. Wind load characteristics of crescent-shaped cantilevered roof with wavy surface[J]. Journal of Zhejiang University(Engineering Science), 2013, 47(6):1063-1071.

[10] DAVENPORT A G. Gust loading factors[J]. Journal of Structural Division, 1967, 93(ST3): 11-34.

[11] KAREEM A, ZHAO J. Analysis of non-Gaussian surge response of tension leg platforms under wind loads[J]. Journal of Offshore Mechanics and Arctic Engineering, 1994, 146: 137-144.

[12] SADEK F, SIMIU E. Peak non-Gaussian wind effects for database-assisted low-rise building design[J]. Journal of Engineering Mechanics, 2002,128: 530-539.

[13] KE S T, JUN L, ZHAO L, et al. Influence of ventilation rate on the aerodynamic interference between two extra-large indirect dry cooling towers by CFD[J]. Wind & Structures an International Journal, 2015, 20(3):449-468.

[14] NISHIMURA H, TANIIKE Y. Aerodynamic characteristics of fluctuating forces on a circular cylinder[J]. Journal of Wind Engineering & Industrial Aerodynamics, 2001, 89(1): 713-723.

[15] CHENG X X, ZHAO L, GE Y J, et al. Wind pressures on a large cooling tower[J]. Advances in Structural Engineering, 2015, 18:201-220.

[16] WINTERSTEIN S R. Nonlinear vibration models for extremes and fatigue[J]. Journal of Engineering Mechanics, 1988, 114(10): 1772-1790.