非穩態雷暴沖擊風場的瞬態數值模擬

汪之松， 王 超， 劉亞南， 李正良

(1.重慶大學 土木工程學院，重慶 400045; 2.重慶大學 山地城鎮建設與新技術教育部重點實驗室，重慶 400045;3.湖南省建筑設計院，長沙 410001)

非穩態雷暴沖擊風場的瞬態數值模擬

汪之松1,2， 王 超3， 劉亞南1， 李正良1,2

(1.重慶大學 土木工程學院，重慶 400045; 2.重慶大學 山地城鎮建設與新技術教育部重點實驗室，重慶 400045;3.湖南省建筑設計院，長沙 410001)

沖擊射流是雷暴沖擊風場研究中最常用的模擬方法。大部分物理試驗和數值模擬都是假定射流速度入口風速不隨時間變化，而實際雷暴沖擊風的下沉氣流速度都是隨時間連續變化的。在雷暴沖擊風的全生命周期內一般會經歷逐漸增大到最大值后再逐漸減小的過程。對于重要的工程結構抗風設計而言，得到雷暴沖擊風全生命周期的風速時程信息十分有必要。基于沖擊射流三維足尺模型模擬雷暴沖擊風風場，在入口處引入一個更符合雷暴沖擊風實際演化過程的衰減函數，采用大渦模擬數值分析獲得了雷暴沖擊風的非穩態風場，并得到其隨時間衰減的瞬態演化過程。結果表明，該模擬方法可以較好地再現雷暴沖擊風場的非穩態過程，為進一步討論非穩態雷暴沖擊風場中的結構風荷載特性奠定了基礎。

雷暴沖擊風；非穩態；沖擊射流；大渦模擬；數值模擬

雷暴沖擊風下沉氣流沖擊地面后，其水平風速在近地面區域產生劇烈的低空風切變，對建筑物和輸電塔等結構有強烈的破壞作用。近年在國內發生了大量的輸電線塔風致倒塌事故，有記錄的75起輸電線塔災害事故中絕大多數都是由雷暴天氣產生的雷暴沖擊風造成[1-2]。國內外眾多學者已采用實際觀測、理論研究、物理實驗和數值模擬等方法對雷暴沖擊風進行了研究，并有一些學者提出了沖擊風徑向和豎向風剖面的解析、經驗模型[3-6]，研究重點從最早的集中于實測記錄和數據收集，到后來的穩態和瞬態的沖擊風理論研究和實驗、數值模擬。SHEHATA[7]研究了一個跨度為480 m的輸電線路沖擊風荷載，得出沖擊風的整個生命周期和強度衰減是非常重要的參數。ABD-ELAAL[8]基于數值模擬研對非穩態雷暴沖擊風的影響因素進行參數分析，表明獲取雷暴沖擊風的非穩態特性及其整個生命周期和強度衰減的信息十分必要。

關于沖擊風的非穩態研究，FUJITA[9]記載了發生于1983年7月1日美國華盛頓圣安德魯斯空軍基地(AAFB)的一次強烈的沖擊風風速時程表現出強烈的非平穩特性。HJELMFELT[10]總結了(Joint Airport Weather Studies,JAWS)研究的數次沖擊風記錄與出流強度隨時間變化情況，表明沖擊風從開始到最大強度幾乎是呈線性增加的，然后再迅速衰減或者在一段時間內強度保持恒定。HOLMES等[11]通過一個衰減函數乘以徑向速度，第一個提出了描述雷暴沖擊風的時程曲線；CHAY等[12-13]更進一步發展了之前的概念，并且引入了幾個強度衰減函數來描述風速時程曲線。CHAY等[14]利用HJELMFELT等[15]的觀測結果給出了一個分析模型的強度衰減函數，他們認為強度在前5～9 min內線性增加，然后在接下來的5～9 min內呈指數衰減(如圖1)。李錦華等[16]對雷暴沖擊風非平穩脈動風速進行了數值合成，引入非均勻調制函數實現非平穩脈動風速的時變功率譜。

本文基于沖擊射流三維足尺模型模擬雷暴沖擊風場，在入口處引入一個更符合雷暴沖擊風實際演化過程的衰減函數，采用大渦(LES)數值模擬獲得了雷暴沖擊風的非穩態風場，并得到其隨時間衰減的瞬態演化過程。

圖1 雷暴沖擊風出流強度隨時間變化的實測記錄[15]Fig.1 Normalized radial velocity versus time[15]

1 CFD數值模擬

用計算流體力學軟件Fluent14.0來模擬不可壓縮的雷暴沖擊風場，采用沖擊射流模型作為雷暴沖擊風的基本風場模型。為了盡可能反映真實物理尺度的雷暴沖擊風現象，建立三維足尺計算域來模擬非穩態沖擊射流場，模型參數包括：初始出流直徑Djet=600 m，計算入口到地面的距離Hjet=2Djet，邊界到噴口中心的徑向距離R=6Djet，計算域示意圖如圖2(a)，(b)所示。

圖2 計算域示意圖Fig.2 The view of the computational domain

為了得到精確的數值模擬結果，采用結構化網格劃分計算域，在沖擊射流中心區域采用4層O型網格，在近壁面區域滿足無量綱距離y+≈1，y+定義為[17]：

(1)

式中:Δy為首層網格到壁面的最小距離，單位m；ν指空氣的運動黏性系數，單位m2/s；τw為壁面切應力，單位Pa；ρ為空氣的密度，單位kg/m3。

計算域網格數量約為1.2×107；速度入口取湍流強度I=1%；頂部及側邊均為壓力出口，壓力為1.013×105Pa；近地面壁面采用增強壁面處理模型；滑移壁面剪應力取為0。

選用大渦模型(LES)對沖擊射流計算域進行求解。該方法將湍流的大渦和小渦分開處理，大漩渦直接求解Navier-Stokes方程，而小渦采用湍流模型來解決。沖擊射流模擬計算的時間步長為0.1 s，總迭代步數6 000步，總計算時長T=10 min。空氣密度為1.25 kg/m3，黏性系數為1.789 4×10-5Ns/m2。

2 速度入口函數

實際的雷暴沖擊風其速度入口出流風速隨時間是連續變化的，一般有一個先逐漸增大到最大值然后再逐漸減小的過程。Abd-Elaal[14]定義了一個速度入口分段函數為：第一段為恒定值，第二段直接減小到接近于0，來描述沖擊風的衰減過程。本文為更加真實地重現雷暴沖擊風的衰減過程，定義依然采用兩段的分段函數來表示速度入口，具體為：第一個時間段Δt1內仍然為恒定值30 m/s；第二個時間段Δt2內為呈指數逐漸減小到0，函數表達式為：

(2)

為了系統地考慮沖擊風衰減周期對整個風場的影響，這里第一個時間段Δt1選取2 min，3 min和4 min三種情況，對應的第二個時間段Δt2=T-Δt1分別為8 min，7 min和6 min，速度入口速度變化函數曲線如圖3所示。

圖3 入口速度變化曲線Fig.3 Inlet velocityversus time

3 非穩態計算結果及分析

3.1 穩定階段平均風速剖面

上述三種情況CFD數值模擬結果的穩定階段(即Δt1段的后30 s風速)水平向平均風剖面與各個分析模型，經驗模型和實測數據豎直風剖面和徑向風剖面進行對比如下圖4、5所示。CFD模擬結果選取r=1.5Djet處的水平風豎直剖面和z=0.01Djet高度處的水平風徑向剖面，與Oseguera和Bowles模型[2]、Vicroy模型[3]、Wood和Kwok模型[4]、Hjelmfelt實測數據[10]和Letchford和Illidge[18]的物理實驗模型等結果進行對比，其中z為豎向高度，zm為最大水平風速對應的豎向高度，u為沿徑向的水平風速，um為沿徑向水平風速最大值，r為相對噴口中心的徑向距離，rm為最大水平風速對應的徑向距離。

圖4 CFD數值模擬與各模型豎直風剖面比較Fig.4 Comparison of vertical wind profile of CFD numerical simulation and analytical models

圖5 CFD數值模擬與各模型徑向風剖面比較Fig.5 Comparison of radial wind profile of CFD numerical simulation and analytical models

從圖中可以看出各個分析模型，經驗模型和實測數據與本文CFD模擬結果都給出了雷暴沖擊風水平風豎直剖面和徑向剖面的主要特征。在水平風豎直剖面方面，由于本文是基于沖擊射流模型的模擬結果，所以和基于沖擊射流模型的物理實驗得出的Wood和Kwok模型基本一致。在水平風徑向剖面方面，CFD數值模擬結果對達到最大水平風速后風速衰減段的風速值過高估計，主要原因是實際的沖擊風經歷形成、持續并最終擴散的過程，并沒有停留在穩定的狀態，CFD數值模擬是提取的風場穩定階段的數據。總體來說，基于LES瞬態風場數值模擬平均水平風剖面同經驗風剖面及實測數據吻合較好，在缺乏實測數據的情況下，除了風洞試驗可以提供較近似的結果外，基于LES的雷暴沖擊風場數值模擬有利于了解其發展全過程尤其是達到最大水平風速前流場的特性。

3.2 流場風速云圖和矢量場

選取Δt1=3 min時的數值分析結果，沖擊射流下沉，發展，穩定和衰減階段的速度云圖和矢量圖如圖6(a)～(f)所示。從圖中可以看到，在t=1 min時，速度入口下沉氣流與周圍氣體之間拖拽卷吸，在下沉氣流周圍形成環形漩渦(a)；在t=2 min時，入口速度保持30 m/s不變，環形漩渦沖擊地面后沿徑向水平發展，風場迅速向周圍輻射擴散，在近地面達到最大出流強度(b)；在3 min過后，入口速度開始減小，到3.5 min時，入口速度下降到最初的50%，整個風場從入口位置開始逐漸衰減，近地面出流強度逐漸減小(c),(d)；從4.5 min到5.5 min，入口速度迅速減小到接近為零，此時整個風場都迅速衰減，隨著近地面出流強度的進一步減弱，雷暴沖擊風逐漸消退(e),(f)。

圖6 沖擊射流各時刻風度云圖和矢量場Fig.6 The contour and vector plot of impinging jet at different time

3.3 水平風速的瞬態分布

選取Δt1=3 min時的數值分析結果來考察雷暴沖擊風整個過程水平風速的瞬態分布情況。

3.3.1 豎直分布

雷暴沖擊風水平風速的豎直分布表征了水平風速沿高度的變化規律，對結構風工程有重要參考價值。為了細致的考察沖擊風整個過程，圖7給出了沖擊射流發展過程中不同徑向位置上水平風速的豎直分布情況。從圖7中可以看到在雷暴沖擊風沖擊地面前后和衰減過程中不同徑向位置水平風速的豎直分布呈現較大差異。隨著下沉氣流沖擊地面并沿地面像四周擴展，水平環形渦流一次經過距離沖擊中心由近至遠的各個徑向位置，導致超過一定高度后出現負向的速度。總體而言，各個徑向位置的最大水平風速在2～3 min左右迅速增加到最大值，然后隨著時間增加逐漸減小，到5 min左右以后在整個豎向位置都下降到最大風速的40%以內。

3.3.2 徑向分布

雷暴沖擊風水平風速的徑向分布表征了水平風速沿徑向的變化規律，圖8給出了沖擊射流發展過程中不同高度位置上水平風速徑向分布情況。從圖8中可以看到，雷暴沖擊風沖擊地面后在結構風工程關心的高度范圍內水平風速的徑向分布具有相似的規律，最大水平風速出現在2 min左右，位于r=1.5Djet位置附近。然后隨著時間的增加，水平風速逐漸減小，到5 min左右以后在整個徑向位置都下降到最大風速的40%以下，這點跟豎直分布所觀察到的一樣，說明雷暴沖擊風在最大風速衰減2 min左右整個風場水平風速都迅速減小到最大值的40%以下。

3.4 豎向風速的瞬態分布

選取Δt1=3 min時的數值分析結果來考察雷暴沖擊風整個過程豎向風速的瞬態分布情況。

3.4.1 豎直分布

雷暴沖擊風豎向風速的豎直分布表征了豎向風速沿高度的變化規律，為了細致的考察沖擊風整個過程，圖9給出了沖擊射流發展過程中不同徑向位置上豎向風速的豎直分布情況。從圖9中可以看到在雷暴沖擊風沖擊地面前后和衰減過程中不同徑向位置豎向風速的豎直分布呈現較大差異。隨著下沉氣流沖擊地面并沿地面像四周擴展，水平環形渦流一次經過距離沖擊中心由近至遠的各個徑向位置，在水平環形渦流經過的徑向位置，無論是風速值還是風向，豎向風速的豎直分布呈現出較大的變化。總體而言，靠近噴口中心的r=0.6Djet在1 min內就很快達到了最大值；r=1.0Djet和r=1.5Djet位置都是在2 min內達到最大值，然后隨著時間增加而逐漸減小；而離噴口中心較遠的r=2.0Djet位置各個時刻的豎向風速都比較小，基本都在0.2Vjet范圍以內，由此可知，豎向風速的影響范圍主要到達1.5Djet左右。

圖7 沖擊射流全過程各徑向位置水平風速的豎直分布Fig.7 Vertical profiles of radial velocity at different time and different radial distances

圖8 沖擊射流過程中各高度位置水平風速的徑向分布Fig.8 Radial profiles of radial velocity at different time and different heights

3.4.2 徑向分布

雷暴沖擊風豎向風速的徑向分布表征了豎向風速沿徑向的變化規律，圖10給出了沖擊射流發展過程中不同高度位置上豎向風速徑向分布情況。從圖10中可以看到，雷暴沖擊風沖擊地面后在結構風工程關心的高度范圍內豎向風速的徑向分布形狀上具有相似的規律，最大豎向風速出現在2 min左右，位于r=1.5Djet位置附近，但最大值的峰值大小有較大差異。然后隨著時間的增加，豎向風速逐漸減小，到5 min左右在整個徑向位置的風速值都下降到0.1Vjet以下。

3.5 水平速度和豎向速度時變特性

數值模擬結果水平速度在Δt1=2,3,4 min三種情況下的時程曲線如圖11所示(徑向位置r=1.0Djet，高度位置z=0.03Djet處)。可以看到三種情況下都有環形渦流經過，導致在2 min左右達到了出流的最大強度，隨著Δt1的增加出流強度沒有更進一步的增加。當入口速度開始下降的時候，監測點速度并沒有突然的下降，而是保持了一段時間的穩定值然后再開始下降進入衰減階段。

圖9 沖擊射流全過程各徑向位置豎向風速的豎直分布Fig.9 Vertical profiles of vertical velocity at different time and different radial distances

圖10 沖擊射流過程中各高度位置豎向風速的徑向分布Fig.10 Radial profiles of vertical velocity at different time and different heights

選取Δt1=3 min時的數值分析結果來考察雷暴沖擊風整個過程風速的時變特性。在徑向r=1.0Djet處結構工程比較關心的高度范圍內(如z=0.002、0.005、0.025、0.1Djet)的水平風速和豎向風速時程曲線如圖12、13所示；在高度z=0.01Djet處各個徑向位置(如r=1.0、1.5、2.0、2.5、3.0Djet)的水平風速和豎向風速時程曲線如圖14和15所示。

圖12展示了速度強度時程曲線隨高度的改變，速度強度剛開始隨著高度的增加而增大，直到最大速度所在高度位置然后減小，其變化趨勢和穩態模型水平風速的豎直風剖面的變化趨勢一致；圖13展示了豎向風速時程曲線隨高度的改變，速度強度隨高度連續的增加，其變化趨勢和穩態模型豎向風速的豎直風剖面的變化趨勢也一致；圖14、15展示了速度強度時程曲線隨徑向距離的改變，速度強度隨徑向距離的增加而增大，到r=1.5Djet處達到最大值然后減小，其變化趨勢與穩態模型水平風速的徑向風剖面也基本一致。

圖11 數值模擬時程曲線Fig．11Temporalprofilesoftheradialspeedforthreedifferentages圖12 r=1．0Djet各高度位置水平風速時程曲線Fig．12Temporalradialspeedatr=1．0Djetfordifferentheights圖13 r=1．0Djet各高度位置豎向風速時程曲線Fig．13Temporalverticalspeedatr=1．0Djetfordifferentheights

圖14 z=0.01Djet處各徑向位置水平風速時程曲線Fig.14 Temporal radial speed at z=0.01Djetfor different radical locations

圖15 z=0.01Djet處各徑向位置豎向風速時程曲線Fig.15 Temporal vertical speed at z=0.01Djet for different radical locations

4 結 論

利用大渦數值模擬方法，基于沖擊射流模型完成了三維足尺的雷暴沖擊風非穩態模擬。在CFD計算結果與相關經驗、分析模型和實測數據進行比較以驗證CFD數值模型及參數準確性的基礎上，進行了風場特性的全面研究。在速度入口處引入一個更符合雷暴沖擊風演化過程的速度衰減函數，得到了沖擊風場整個生命周期及衰減過程的詳細信息。結構風工程關心的沖擊風近地面風場的水平風速和豎向風速總體而言都在短時間內迅速達到最大值，隨后逐漸衰減或保持一個恒定值波動后再迅速衰減至零。基于沖擊射流模型開展的CFD瞬態數值模擬可以較好地體現雷暴沖擊風整個衰減過程的基本特性，為進一步討論非穩態雷暴沖擊風場中的結構風荷載奠定了基礎。

[ 1 ] 張勇. 輸電線路風災防御的現狀與對策[J]. 華東電力, 2006,34(3):28-31. ZHANG Yong. Status quo of wind hazard prevention for transmission lines and countermeasures [J]. East China Electric Power, 2006,34(3):28-31.

[ 2 ] 瞿偉廉，吉柏鋒. 下擊暴流的形成與擴散及其對輸電線塔的災害作用[M]. 北京：科學出版社，2013.

[ 3 ] OSEGUERA R M, BOWLES R L. A simple analytic 3-dimensional downburst model based on boundary layer stagnation inflow [R].NASA Technical Memorandum 100632, 1988.

[ 4 ] VICROY D D.A simple, analytical, asymmetric microburst model for downdraft estimation [R]. NASA Technical Memorandum 104053, 1991.

[ 5 ] WOOD G S, KWOK K C S, An empirically derived estimate for the mean velocity profile of a thunderstorm downbursts [C]//Proceedings of 7th AWES Workshop. Auckland, 1988.

[ 6 ] HOLMES J D, OLIVER S E. An empirical model of a downburst [J].Engineering Structures, 2000, 22: 1167-1172.

[ 7 ] SHEHATA A Y, EL DAMATTY A A, Savory E. Finite element modeling of transmission line under downburst wind loading [J]. Finite Elements in Analysis and Design, 2005, 42(1):71-89.

[ 8 ] ABD-ELAAL E, MILLS J E, MA X. A coupled parametric-CFD study for determining ages of downbursts through investigation of different field parameters [J]. Journal of Wind Engineering and Industrial Aerodynamics, 2013, 123: 30-42.

[ 9 ] FUJITA T T. The downburst report of Projects NIMROD and JAWS [R]. Chicago: University of Chicago, 1985.

[10] HJELMFELT M R. Structure and life circle of microburst outinflows observed in Colorado[J]. J. ApplMeteorol, 1988, 27(8): 900-927.

[11] HOLMES J D, OLIVER S E. An Emprical Model of a Downburst[J].Engineering Structures, 2000, 22:1167-1172.

[12] CHAY M T, ALBERMANI F, WILSON R. Numerical and analytical simulation of downburst wind loads [J]. Engineering Structures, 2006, 28(2): 240-254.

[13] ABD-ELAAL E, MA X, MILLS J E. A newly developed analytical model of transient downburst wind loads [C]//Proceedingsof the 22th Australian Conference on the Mechanics of Structures and Material, Sydney, Australia, 2012: 425-430.

[14] ABD-ELAAL E, MILLS J E, MA X. Empirical models for predicting unsteady-state downburst wind speeds [J]. J. Wind Eng. Ind. Aerodyn, 2014, 129: 49-63.

[15] WILSON J W, ROBERTS R D, KESSINGER C, et al. Microburst wind structure and evaluation of Doppler radar for airport wind shear detection [J]. Journal of Climate and Applied Meteorology, 1984, 23(6):898-915.

[16] 李錦華，吳春鵬，陳水生.下擊暴流非平穩脈動風速數值模擬[J]. 振動與沖擊，2014, 33(14):54-60. LI Jinhua，WU Chunpeng，CHEN Shuisheng. Simulation of non-stationary fluctuating wind velocity in downburst[J]. Journal of Vibration and Shock,2014, 33(14):54-60.

[17] Fluent Inc. FLUENT User’s Guide[Z]. Fluent Inc., 2003.

[18] LETCHFORD C W, ILLIDGE G. Turbulence and topographic effects in simulated thunderstorm downdrafts by wind tunnel jet[C]//Proceedings of the Tenth International Conference on Wind Engineering, Denmark, June 1999:1907-1912.

Transient numerical simulation for unsteady-state downburst winds

WANG Zhisong1,2,WANG Chao3, LIU Yanan1, LI Zhengliang1,2

(1. School of Civil Engineering, Chongqing University, Chongqing 400045, China;2. Key Laboratory of New Technology for Construction of Cities in Mountain Area, Ministry of Education, Chongqing University, Chongqing 400045, China;3. Hunan Provincial Architectural Design Institute, Changsha 410001, China)

Impinging jet is one of the most commonly used simulation methods in studying the wind field of thunderstorms. Most of the previous physical tests and numerical simulations assume that the inlet velocity of the jet does not change with time, while the sinking velocity of downburst winds is continuously changing with time actually during its whole lifecycle, it gradually increases to the maximum value and then decreases in general. For wind-resistant design of significant engineering structures, it’s necessary to acquire the temporal profiles of downburst wind speed in its whole lifecycle. Here, a three-dimensional full-scale impinging jet was taken as a wind field model of downburst. The decay functions for the inlet jet velocity were given for the actual downburst，and numerical simulation was performed by adopting the large eddy simulation(LES) to acquire an unsteady-state wind field where the intensity decay information of downburst was presented. The simulation results showed that the unsteady wind field of downburst can be well simulated with this proposed method, and the results are correlated well to Andrews AFB recorded ones. These results laid a foundation for further studying wind load characteristics of unsteady downburst.

downburst; unsteady-state; impinging jet; LES; numerical simulation

國家自然科學基金(51208537)

2015-10-09 修改稿收到日期：2016-01-27

汪之松 男，副教授，1980年生

TH212；TH213.3

A

10.13465/j.cnki.jvs.2017.03.009