祝志文，袁濤，陳魏

(湖南大學 土木工程學院，湖南 長沙 410082)

?

扁平箱梁氣動特性CFD模擬的維數(shù)對比研究

祝志文，袁濤，陳魏

(湖南大學 土木工程學院，湖南 長沙 410082)

基于三維非定常雷諾時均Navier-Stokes (RANS)方程和SST k-ω湍流模型，數(shù)值模擬雷諾數(shù)為1.06×105時扁平箱梁氣動特性在不同主梁展向長度、不同展向網(wǎng)格分辨率下的氣動特性，確定三維CFD模擬合理的主梁展向參數(shù)，獲得主梁氣動力系數(shù)和漩渦脫落St數(shù)隨來流攻角的變化。對相同箱梁開展了不同攻角下二維RANS計算，并將二維和三維主梁的整體氣動特性、表面壓力分布進行對比。研究結(jié)果表明：三維RANS獲得的主梁氣動特性結(jié)果與計算域展向網(wǎng)格分辨率和主梁長度均無關(guān)；二維和三維模擬結(jié)果非常一致，均能獲得箱梁氣動力系數(shù)平均值的合理估計，因而工程計算可采用二維模擬。開展的節(jié)段模型風洞試驗表明二維和三維模擬得到的漩渦脫落St數(shù)是合理的，典型攻角下扁平箱梁裸梁具有單一St數(shù)，且均明顯大于成橋狀態(tài)。

扁平箱梁；氣動特性； RANS；二維和三維模擬

湍流是一種極其復(fù)雜的流動形式。為獲得對湍流流動的描述，雷諾在1886年就提出了用時間平均值的概念來研究湍流流動。他認為湍流中的所有物理量雖然都隨時間和空間變化，但任意瞬時的運動仍然符合連續(xù)介質(zhì)流動的特征，流場中任一空間點也應(yīng)適用流體運動的基本方程。且由于各個物理量都具有某種統(tǒng)計學特征，所以控制方程中任一瞬時的物理量都可用其平均量和脈動量之和來代替，因而可對整個Navier-Stokes方程進行時間平均運算。

對不可壓流動，任一瞬時的流動物理量可表示為：

(1)

(2)

(3)

其中：t和x分別為時間和空間坐標；ρ和υ分別為空氣密度和運動黏性。

動量方程可進一步寫成

(4)

風工程中幾乎所有流動問題，由于流動雷諾數(shù)高，其繞流流動均是復(fù)雜的湍流流動。因而目前風工程數(shù)值模擬需要借助湍流模型來獲得流場描述、風荷載大小和加載特征。又由于風工程中的繞流物體往往具有鈍化的氣動外形和尖銳的棱角，流動總是出現(xiàn)強烈的分離、物面再附、再分離、漩渦脫落和復(fù)雜的尾流流動。因此所采用的湍流模型應(yīng)能處理高應(yīng)變率及大彎曲程度的流動，能模擬逆壓和分離流動，具有較強的捕捉流動各向異性的能力。本文采用SSTk-ω湍流模型研究橋梁斷面的氣動特性，并依此開展二維和三維模擬的比較。

1　湍流模型方程

如果采用渦黏模型封閉上述方程, 對應(yīng)的雷諾應(yīng)力可表示為,

(5)

Menter提出了一種SST(ShearStressTransport)k-ω模型[2]。該模型通過對SSTk-ω模型方程和標準k-ε模型方程分別乘以混合函數(shù)Fl和(1-Fl)。通過設(shè)定近壁區(qū)內(nèi)Fl取值1，使得近壁區(qū)的流動采用SSTk-ω模型模擬；在近壁區(qū)以外Fl為0，從而過度到標準k-ε模型。其中標準k-ε模型方程為[3]

(6)

(7)

其中:σk和σε分別是k和ε的湍流普朗特數(shù)；Cε1和Cε2為常數(shù)。標準k-ε模型是一種高雷諾數(shù)湍流模型，適合于充分發(fā)展的湍流，但在湍流邊界層內(nèi)緊靠壁面區(qū)的黏性子區(qū)內(nèi)，流動雷諾數(shù)低，并不是充分發(fā)展的湍流，需要使用如下的低雷諾數(shù)k-ω湍流模型[4]：

(8)

(9)

這樣， SST k-ω湍流模型結(jié)合了標準k-ε模型和k-ω湍流模型各自的優(yōu)缺點，通過不同湍流模型用于不同的流動區(qū)域計算，能獲得優(yōu)于標準k-ε模型和SSTk-ω模型的模擬結(jié)果[5]。

2　研究對象

以丹麥大帶東橋主跨加勁梁為研究對象。該橋為主跨1 624 m的雙塔三跨連續(xù)鋼箱梁懸索橋，加勁梁橫斷面為扁平閉口六邊形，橋面布置如圖1所示。該加勁梁全寬B=31 m，橋軸線處梁高4.4 m。該橋在設(shè)計階段開展了大量的風洞試驗[6]。為與風洞試驗對比，CFD模擬采用與節(jié)段模型風洞實驗一致的縮尺比和來流風速[6]，對應(yīng)的流動Re=1.06×105。由于作者無法獲得橋梁欄桿、導(dǎo)流板和檢修道等附屬設(shè)施的詳細結(jié)構(gòu)設(shè)計細節(jié)，本文CFD計算不考慮橋面欄桿和其他附屬設(shè)施，相當于加勁梁的施工狀態(tài)。

圖1　加勁梁橫斷面布置Fig.1 Layout of girder cross section

3　數(shù)值方法實現(xiàn)

二維和三維CFD模擬均采用1:80的模型縮尺比，繞加勁梁橫斷面的計算域如圖2所示。計算域的入口、上側(cè)和下側(cè)邊界到加勁梁斷面中心的距離均為13B，對應(yīng)的模型堵塞度為0.5%，顯著小于風工程對模型堵塞度不大于3%的要求[7]。下游出口到加勁梁扭轉(zhuǎn)中心(Shear center，S.C.)的距離為26B，目的是減小下游邊界對模型區(qū)流動的影響。采用計算域分區(qū)劃分網(wǎng)格以控制網(wǎng)格的正交性和相鄰網(wǎng)格邊長比不大于1.2，計算域分區(qū)網(wǎng)格劃分見文獻[8]。由于前后緣尖角小半徑圓化對主梁氣動特性沒有影響[8]，為減少高Re數(shù)下的計算域網(wǎng)格總數(shù)，提高主梁前后緣網(wǎng)格正交性，前后緣圓化半徑均為0.25 mm，如圖3所示。經(jīng)計算結(jié)果的網(wǎng)格無關(guān)檢查，確定二維主梁物面第一條網(wǎng)格到物面的最大距離為2.9×10-5B，計算Re數(shù)下的最大Y+為0.43，以度量網(wǎng)格正交性的等角斜度(Equiangle skew)QEAS參數(shù)最大值為0.4，保證了較好的網(wǎng)格正交性。二維計算域網(wǎng)格總數(shù)為55 332，其他詳細信息見文獻[8]。對展向不同主梁長度、不同網(wǎng)格分辨率的三維計算模型，其任意主梁橫截面的網(wǎng)格劃分和參數(shù)與二維網(wǎng)格完全相同，如圖4所示。

圖2　加勁梁與計算域Fig.2 Schematic plot of computational domain around girder

計算域邊界條件為：計算域入口為水平層流速度邊界，下游出口施加牛曼邊界條件；加勁梁表面為無滑移壁面條件；計算域上、下邊界采用對稱邊界條件，三維主梁的2個端截面所在邊界同樣采用對稱邊界條件。計算域初始場采用入口速度條件初始化。數(shù)值計算采用循序漸進的推進策略，以平抑初始場計算導(dǎo)致的數(shù)值振蕩和數(shù)值發(fā)散，縮短計算時間。最終的CFD計算控制方程時間離散采用非定常二階隱式格式，采用速度-壓力解耦的SIMPLE算法，壓力方程采用二階格式離散，對流項、湍動能方程和湍流耗散率方程均采用QUICK格式。經(jīng)二維計算的時間步無關(guān)檢查，確定的無量綱時間步長Δ=0.003 1，并用于三維CFD計算。通過監(jiān)視氣動力時程，當加勁梁上作用的氣動力收斂后開始保存每步上的氣動力和壓力數(shù)據(jù)，所有數(shù)值模擬均基于CFD專用程序Fluent 6.3開展。

圖3　加勁梁前緣網(wǎng)格Fig.3 Grids around girder leading edge

圖4　三維計算域主梁橫斷面網(wǎng)格Fig.4 Mesh slice at girder cross section in three-dimensional domain

4　數(shù)值模擬結(jié)果和討論

4.1計算結(jié)果的處理和監(jiān)測

定義加勁梁斷面氣動升力、阻力和扭矩系數(shù)分別為：

(10)

其中:U0為計算域入口風速,FL，F(xiàn)D和M分別為作用在加勁梁上的阻力(順來流方向為正)、升力(垂直來流方向向上為正)和扭矩(逆時針方向為正)；H和L分別為加勁梁斷面高度和模擬的加勁梁展向長度，二維時L為單位長度。

定義加勁梁漩渦脫落St數(shù)：

St=fsH/U0

(11)

式中:fs為加勁梁漩渦脫落頻率,Hz。

為獲得主梁表面的壓力分布特征，本文分別在二維主梁截面和三維主梁展向中截面上相同位置布設(shè)了160個監(jiān)測點。壓力監(jiān)測點的分布疏密根據(jù)表面壓力梯度的大小調(diào)整，在前后緣點位置加密，如圖5所示。定義測點的壓力系數(shù)為

(11)

其中:P為監(jiān)測點靜壓；Pr為參考點壓力，參考點均設(shè)在計算域入口邊界。

4.2二維CFD模擬結(jié)果

二維模擬的來流攻角范圍為-6～+6，得到的氣動力系數(shù)如圖6所示?？梢娫谀M的攻角范圍內(nèi)，阻力系數(shù)0°攻角最小，偏離0°攻角均逐漸增大；升力系數(shù)隨攻角的增大而增大，扭矩系數(shù)隨攻角的增大而減小。另外，0°攻角下的阻力系數(shù)和升力系數(shù)，均與不帶欄桿截面二維模擬采用RNGk-ε模型和雷偌應(yīng)力模型(RMS)接近[9]。風洞試驗阻力系數(shù)明顯偏大，其原因是欄桿和附屬設(shè)施明顯增大了截面阻力[6]。且。圖7是St數(shù)隨來流攻角的變化，其值在0.27～0.28，表現(xiàn)為對來流攻角變化不敏感。圖8是0°攻角下主梁三分力系數(shù)時程，可見計算結(jié)果已充分收斂，力系數(shù)時程已非常平穩(wěn)。雖然時程脈動量非常小，但仍然可基于升力時程的渦脫頻率(圖9)獲得渦脫St數(shù)的估計值。

圖5　主梁斷面壓力監(jiān)測點布置Fig.5 Arrangement of pressure monitoring point on girder surface

圖6　三分力系數(shù)隨來流攻角的變化Fig.6 Force coefficients vs. wind angle of attack

圖7　主梁St數(shù)隨來流攻角的變化Fig.7 Girder St number vs. wind angle of attack

圖8　主梁三分力系數(shù)時程Fig.8 Force coefficients records of girder

圖10是0°攻角下主梁表面平均壓力系數(shù)分布(曲線與物面間畫填充線的為負壓，后同)?？梢姵熬壣舷赂拱迳喜糠謪^(qū)域為正風壓外，主梁表面其它區(qū)域均為負風壓，在頂?shù)装遄钌嫌蝹?cè)的棱角，負壓峰值很大。最大正風壓系數(shù)0.95位于前緣點下，最大負風壓系數(shù)-1.94位于頂板最上游點，下游表面的風壓系數(shù)一般在-0.3～-0.2范圍，主梁大部分區(qū)域風壓系數(shù)絕對值較小。

圖9　升力系數(shù)時程功率譜分析Fig.9 PSD plot of lift coefficient records

圖10　主梁表面平均壓力系數(shù)分布Fig.10 Pressure coefficients distribution on girder surface

4.3三維CFD模擬結(jié)果

對橋梁主梁的三維CFD模擬，需要研究主梁展向網(wǎng)格分辨率和主梁展向長度對其氣動特性的影響，由此確定主梁氣動特性模擬的合理展向長度和網(wǎng)格分辨率。

4.3.1主梁展向網(wǎng)格分辨率對氣動特性計算結(jié)果的影響

針對主梁展向長度為1倍主梁寬度的情況，計算了0° 攻角下5種不同主梁網(wǎng)格分辨率，如表1所示。對不同主梁展向網(wǎng)格分辨率，在垂直主梁軸向的平面，網(wǎng)格劃分與二維模擬完全相同。主梁氣動力系數(shù)和渦脫St數(shù)見表2，同時列出了二維RANS的相應(yīng)參數(shù)值。

表1顯示5種不同網(wǎng)格分辨率的氣動三分力系數(shù)完全相同，而St數(shù)差別也非常小。這說明，采用三維RANS模擬主梁氣動特性與計算域展向網(wǎng)格分辨率無關(guān)。同時，0° 攻角下三維計算結(jié)果與二維結(jié)果非常接近，偏差最大的為阻力系數(shù)，而相對誤差僅為3.9%。

表1不同展向網(wǎng)格分辨率的流動結(jié)果

Table 1 Flow results of girder under various spanwise grid resolution

模型編號展向網(wǎng)格數(shù)量CDCLCMStR-1100.30790.0437-0.0330.281R-2200.30790.0437-0.0330.283R-3300.30790.0498-0.0330.280R-4400.30790.0498-0.0330.277R-5600.30790.0498-0.0330.279二維模擬0.31990.0464-0.0330.279

4.3.2裸梁St數(shù)的節(jié)段模型風洞試驗

為驗證本文CFD模擬獲得的St合理性，在湖南大學HD-2風洞中開展了大帶東橋節(jié)段模型風洞試驗。節(jié)段模型幾何縮尺比為1∶60，模型截面寬B=51.67 cm，高D=7.34 cm，模型總長250 cm，模型為不帶欄桿和附屬設(shè)施的裸梁。模型直立在風洞內(nèi)，在節(jié)段模型的正下游1.5B，離風洞底板0.6 m高度處布設(shè)Cobra-100風速儀，如圖11所示。通過測量模型尾跡縱向脈動風速，由頻譜分析來估算模型的漩渦脫落頻率和St數(shù)。FFT分析點數(shù)為1 024，采用分段平均周期圖法估計功率譜。

圖11　節(jié)段模型試驗Fig.11 Sectional model wind tunnel test

圖12是試驗來流風速分別為5和8 m/s，裸梁在3種典型攻角下尾跡風速測量得到的渦脫St數(shù)，約為0.27～0.29。這說明典型攻角下大帶橋裸梁具有單一St數(shù)，其值為0.28～0.29，均明顯大于成橋狀態(tài)，因此施工階段大跨度橋梁可能會呈現(xiàn)明顯不同于成橋狀態(tài)的渦振特性。

(a)5 m/s;(b)8 m/s圖12　裸梁基于尾跡風速測量的渦脫St數(shù)Fig.12 Vortex shedding St number of bare deck based on wind speed measured in wake

4.3.3主梁展向長度對氣動特性計算結(jié)果的影響

三維CFD模擬需檢查不同主梁展向長度對氣動參數(shù)的影響?；谥髁簹鈩犹匦圆皇苡嬎阌蛘瓜蚓W(wǎng)格分辨率影響的結(jié)論，下面針對計算域展向30個網(wǎng)格的模型，分別模擬計算域5種不同主梁展向長度，得到的氣動參數(shù)見表2。

從表2可見，采用不同的主梁展向長度，其氣動三分力系數(shù)和漩渦脫落St數(shù)幾乎完全一致，且與二維結(jié)果非常接近，最大偏差為阻力系數(shù)，但其相對誤差僅為3.9%?？梢?，三維RANS模擬的氣動特性結(jié)果，不受主梁展向長度的影響。

因三維RANS模擬的主梁氣動特性，與計算域展向網(wǎng)格分辨率和主梁展向長度均無關(guān)，因而三維RANS模擬可采用更為簡化的CFD網(wǎng)格模型。本文確定采用展向1倍主梁寬度和展向10個均勻網(wǎng)格劃分的R-1計算模型，據(jù)此研究在不同來流攻角下的主梁氣動特性。

表2展向不同長度主梁的流動結(jié)果

Table 2 Flow results of girder under various spanwise girder length

模型編號展向長度/BCDCLCMStL-10.50.30790.0497-0.03260.281L-210.30790.0498-0.03260.280L-31.50.30790.0497-0.03260.280L-420.30790.0496-0.03260.282L-530.30790.0497-0.03260.281二維模擬0.31990.0464-0.0330.279

4.3.4不同攻角下主梁氣動特性

三維計算的來流攻角范圍為-6～+6，得到的氣動力系數(shù)和漩渦脫落St數(shù)見圖6～7?？梢娙S結(jié)果與二維結(jié)果具有完全相同的趨勢，且同一攻角下的結(jié)果差別非常小。最大偏差為阻力系數(shù)，但相對差別僅為4%左右。St數(shù)取值同樣在0.28～0.29左右，且表現(xiàn)為對來流攻角變化不敏感。圖13是0°攻角下三維主梁軸向中截面表面平均壓力系數(shù)分布。與圖10相比，平均壓力分布的特征是完全相似的，頂?shù)装遄钌嫌蝹?cè)的棱角負壓峰值均很大，但最小負壓系數(shù)值-2.11小于二維情況，且已改變到位于下底板最上游棱角處。圖14是6°攻角下三維主梁的三分力系數(shù)時程，可見力系數(shù)計算結(jié)果已收斂，但與圖7相比由于攻角增大，力系數(shù)平均值均已增大，且波動幅度也增大。

風洞試驗給出了來流風速和縱向湍流度分別為12 m/s和7.5%，該加勁梁成橋狀態(tài)0°攻角的阻力、升力和扭矩系數(shù)分別為0.57，0.067和0.028[6]。由于本文來流為均勻流，且風洞試驗?zāi)Ｐ桶瑯蛎鏅跅U、防撞墻和檢修道等橋面附屬設(shè)施，因而本文結(jié)果與風洞試驗結(jié)果沒有嚴格的可比性。當防撞欄等橋面附屬設(shè)施引入后，主梁氣動外形將明顯鈍化，小攻角時其阻力系數(shù)可能會增大60%～100%左右，這樣本文二維和三維的阻力系數(shù)值將比較接近風洞試驗結(jié)果。另外，該橋?qū)崢虬l(fā)生渦振的St數(shù)約為0.1～0.14，顯著低于本文計算結(jié)果。其中原因尚不清，可能的原因包括防撞欄等橋面附屬設(shè)施引入后，主梁氣動外形將明顯鈍化，其渦脫特征可能明顯改變；且由于RANS模型只能獲得非定常流動的時均結(jié)果，而漩渦脫落是橋梁繞流流動分離、再附、在分離和剪切流融合的強烈非定常運動，二維和三維RANS模型可能無法獲得漩渦脫落St數(shù)的合理估計。

圖13　主梁中截面平均壓力系數(shù)分布Fig.13 Pressure coefficients distribution on girder span central section

圖14　6°攻角主梁力系數(shù)時程Fig.14 Force coefficients records of girder under 6° angle of attack

圖15為來流攻角為6°時，繞三維箱梁中截面的流線，可見繞箱梁的流動比較平順，在箱梁后緣存在由于流動分流形成的回流區(qū)和上下表面分流剪切流的融合，這些特征同樣可在來流攻角為6°時的二維模擬中見到，如圖16所示。

圖15　6°攻角繞三維箱梁中截面的流線Fig.15 Streamlines around three-dimensional girder at span central section under 6° angle of attack

圖16　6°攻角繞二維箱梁的流線Fig.16 Streamlines around two-dimensional girder under 6° angle of attack

5　結(jié)論

1)基于二維和三維模擬RANS模擬，獲得箱梁氣動力系數(shù)平均值的合理估計。

2)節(jié)段模型風洞試驗和二維、三維模擬得到的漩渦脫落St數(shù)是一致的，典型攻角下扁平箱梁裸梁具有單一St數(shù)，且均明顯大于成橋狀態(tài)。

3)三維RANS計算的主梁氣動特性結(jié)果，與計算域展向網(wǎng)格分辨率和主梁展向長度均無關(guān)，因而三維RANS模擬可采用更為簡化的CFD模型。

4)二維和三維模擬獲得的箱梁氣動特性差別很小，從工程應(yīng)用的精度要求和研究周期，扁平閉口箱梁可僅開展二維CFD模擬，有三維外形特征的主梁除外(如含展向變化的橋面附屬設(shè)施)。

[1] Tannehill J C, Anderson D A, Pletcher R H. Computational fluid mechanics and heat transfer[M]. 2nd edition. Taylor & Francis, Washington, D.C., 1997: 165-176.

[2] Menter F R. Two-equation eddy-viscosity turbulence models for engineering applications[J]. AIAA Journal, 1994, 32: 269-289.

[3] Jones W P, Launder B E. The prediction of laminarization with a two-equation model of turbulence[J]. International Journal of Heat and Mass Transfer, 1972, 15: 301-314.

[4] Wilcox D C. Re-assessment of the scale-determining equation for advanced turbulence models[J]. AIAA Journal, 1988, 26: 1414-1421.

[5] 祝志文,夏昌.基于兩種湍流模型的橋梁顫振導(dǎo)數(shù)識別研究及比較[J].湖南大學學報(自然科學版),2010, 37(11)：6-11.

ZHU Zhiwen, XIA Chang. Comparative study of two turbulence models based on the identification of flutter Derivatives of Bridges[J]. Journal of Hunan University (Natural Sciences), 2010, 37(11): 6-11.

[6] Larsen A. Aerodynamic aspects of the final design of the 1624m suspension bridge across the Great Belt[J]. J Wind Eng Ind Aerodyn, 1993, 48: 261-285.

[7] Tominagaa Y, Mochidab A, Yoshiec R, et al. AIJ guidelines for practical applications of CFD to pedestrian wind environment around buildings[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2008, 96: 1749-1761.

[8] 祝志文.均勻流中箱梁非定常風壓的頻域特征和空間相關(guān)性[J].中國公路學報，2014,27(6):68-76.

ZHU Zhiwen. Frequency-domain characteristics and spatial coherence of unsteady wind pressure on box girder surface in smooth flow[J]. China Journal of Highway and Transport，2014,27(6):68-76.

[9] Bruno L. The validity of 2D numerical simulations of vortical structures around a bridge deck[J]. Mathematical and Computer Modelling, 2003,37(3):795-828.

Comparative investigation on prediction of aerodynamic characteristics of steel flat box girders by CFD with different simulation dimension

ZHU Zhiwen, YUAN Tao, CHEN Wei

(School of Civil Engineering, Hunan University, Changsha 410082, China)

Based on three-dimensional unsteady Reynolds averaged Navier-Stokes (RANS) equations and the SST k-ω turbulent model, the aerodynamic characteristics of a flat box girders under various spanwise length and grid resolution were numerically computed at the Reynolds number of 1.06×105. According to the determined optimal spanwise parameters in CFD, the aerodynamic force coefficients and vortex-shedding Strouhal number were obtained under different wind angles of attack. Flow results of the girder under the same condition were also simulated on the basis of the two-dimensional RANS equations, its aerodynamic forces and pressure on girder surface were compared to the three-dimensional modeling. It is found that the results obtained from three-dimensional simulations is independent on grid resolution and girder length in spanwise direction, and agree well with the results from two-dimensional simulations. Both the two- and three-dimensional RANS modeling can present reasonable estimation on mean aerodynamic coefficients of the girder. Therefore, the two-dimensional modeling is feasible for bridge engineering application. The wind tunnel tests provided confirmed that the estimated vortex-shedding Strouhal number from two- and three-dimensional RANS modeling, and under typical wind angle of attack(-3°～+3°), the bare deck presents only one Stnumber, which is significantly larger than that of girder in service.Key words：box girder; aerodynamics; RANS; two- and three-dimensional modelling

2015-10-14

國家重點基礎(chǔ)研究發(fā)展計劃(973計劃)項目(2015CB057701，2015CB057702)；國家自然科學基金資助項目(51278191)；湖南省交通科技資助項目(201522)

祝志文(1968—)，男，教授，博士，從事大跨度橋梁和結(jié)構(gòu)抗風、鋼橋疲勞與斷裂、橋梁沖刷等研究；E-mail: zwzhu@hnu.edu.cn

TU328

A

1672-7029(2016)08-1555-08