大跨索承橋梁流線型鋼箱梁抖振響應流固耦合數值模擬

徐梓棟 王浩 劉震卿

摘要 大跨索承橋梁抖振性能評估是橋梁抗風安全設計的重要環節。但當前抖振分析理論在揭示抖振基本物理成因時的“描述性”強于“解釋性”，且風洞試驗難以重現抖振過程中的流固耦合細節。對此，以蘇通大橋標準主梁節段為研究對象，在Fluent軟件中開展了模型抖振響應數值模擬，并將抖振響應數值解與理論解進行對比。結果表明，對高湍流度風場而言，模型風振以湍流引起的強迫振動為主，考慮流固耦合效應的抖振響應理論解與數值解的RMS較為接近，模型上、下方渦核發展交替進行，特征湍流對模型抖振響應的影響較小。

關鍵詞 大跨橋梁; 抖振響應; 計算流體力學; 流固耦合; 大渦模擬

引 言

大跨索承橋梁抖振是由空氣中湍流引起的限幅強迫振動，對大跨橋梁抖振性能進行評估是橋梁抗風安全設計的重要環節［1?2］。目前，工程中常用的抖振分析理論均引入了大量的氣動參數對抖振荷載進行建模，且各類氣動參數需由風洞試驗或數值模擬等方法確定，因而既有抖振分析理論被認為具有半試驗性質與統計性質［3?4］。此外，既有抖振分析理論以最大程度刻畫橋梁抖振現象并服務于工程實際為目的，在工程領域得到廣泛應用，但其難以揭示抖振現象的基本物理成因，且實際中重要的流固耦合（Fluid?Structure Interaction， FSI）細節容易被忽視［3，5?6］。因而，有必要發展更加精細化的大跨橋梁抖振響應分析方法，盡量避免使用試驗氣動參數與經驗模型。

為充分把握空氣與橋梁結構的FSI細節，實現大跨橋梁氣動響應的精細化模擬，研究者們通過結合計算流體動力學（Computational Fluid Dynamics， CFD）與結構動力學方法，對橋梁?空氣耦合系統動力方程進行求解［6］。例如，Nariman［7］采用雙向FSI對大跨索承橋梁的渦激振及其“鎖定”現象進行了分析；詹昊等［8］基于用戶自定義函數（User Defined Function， UDF）對Fluent軟件進行二次開發，實現了大跨橋梁主梁顫振的FSI數值仿真，并據此開展了橋梁顫振主動控制研究；董國朝等［9］依托ANSYS CFX軟件及User Fortran工具開發了流固耦合分析程序，并據此開展了斜拉橋拉索風致馳振FSI分析；Kavrakov等［10］采用二維模型開展了大跨橋梁主梁節段抖振響應FSI分析，并將數值解與理論解進行了對比。總體上，現有大跨索承橋梁數值模擬研究多針對如顫振、渦振等自激振動，而受大氣邊界層湍流難以準確模擬、數值計算開銷巨大、氣固耦合機理難以準確把握等因素制約［10?11］，對由湍流引發的大跨橋梁抖振響應數值模擬研究并不多見。

大跨索承橋梁風致抖振是由大氣邊界層湍流引起的強迫振動，因而合理高效地模擬大氣邊界層湍流風環境是開展風致抖振分析的重要前提。大渦模擬（Large Eddy Simulation， LES）是CFD中重要的湍流處理方式，其對大渦部分進行直接求解，而對亞尺度渦旋進行模化。利用LES開展大氣邊界層湍流風環境模擬的關鍵在于生成滿足特定條件的入口湍流［11?12］。當前，LES入口湍流生成技術包括層?湍流過渡技術、循環域技術、驅動域技術與合成湍流技術四類［13］。其中，合成湍流技術又包括譜方法、渦方法與濾波法三類，因計算開銷小且可實現入口邊界平均風剖面與湍流特性的精準控制，該技術被廣泛應用于湍流大氣邊界層模擬［14］。綜合考慮各類LES湍流入口生成技術在計算效率、湍流結構可控性、大氣邊界層模擬適用性等方面的優缺點，擬采用合成湍流技術中的窄帶合成（Narrowband Synthetic Random Flow Generation， NSRFG）法［14］生成滿足給定風特性的湍流風場，并據此實現主梁節段抖振響應的數值模擬。

本文依托Fluent軟件平臺，基于Fluent中的用戶自定義函數（User Defined Function， UDF）功能植入NSRFG法，從而生成滿足給定風特性的湍流風場。此外，通過UDF編制蘇通大橋主梁標準節段的抖振響應FSI模擬程序，開展該橋三維主梁節段抖振響應數值模擬。為確保模擬結果合理可靠，以大橋鋼箱梁氣動參數的風洞試驗結果為參考，對數值模型進行驗證，并將抖振響應數值解與半經驗模型獲得的理論解進行對比分析。

1 蘇通大橋標準主梁節段數值模型

1.1 幾何建模與網格劃分

取模型幾何縮尺比1∶50，并參考文獻［15］構造三維計算域，其中計算域z方向厚度16 m為蘇通大橋標準主梁節段長度，計算域x?y平面尺寸及邊界條件設置如圖1所示。

圖1中，鋼箱梁重心位置坐標為（x， y=0 mm， z=160 mm）。為便于抖振響應數值模擬，計算域網格劃分采用多域網格劃分方案以降低對模擬時間步長的要求。各網格區域由interior交界面分隔，各區域網格在交界面上共節點。在充分利用計算資源的前提下，計算域采用幾何適應性較強的非結構網格劃分，邊界層采用三棱柱網格劃分，第一層邊界層網格厚度為0.001 m，壁面y+平均值45，邊界層中的黏性亞層與緩沖層則通過壁面函數進行模化，計算域網格總數超460萬（4611692個），動網格區域及邊界層網格劃分情況如圖2所示。

圖2中網格劃分在Hypermesh軟件中完成，計算域網格最大等效扭曲度（Equivalent Skew）小于0.7，滿足Fluent計算要求。

1.2 數值模型氣動參數分析與驗證

節段模型的風致抖振主要由抖振荷載引起，抖振理論表明抖振荷載與湍流風場及斷面三分力系數密切相關［10］。其中，湍流風場可通過LES實現，數值模型的三分力系數則需要依據風洞試驗結果進行驗證。

數值模型三分力系數分析時，關閉動網格設置，采用雷諾平均N?S方程（Reynolds?Averaged Navier?Stokes， RANS）進行瞬態（Unsteady）計算，湍流模型取Realizable k ?ε模型，設置變量與連續性方程殘差小于10–5為收斂標準，時間步長取0.002 s，時間步總數13000步，離散與插值格式設置如表1所示。以1.1節中數值模型為分析對象，分別對模型-3°到+3°攻角下（步長1°）的三分力系數進行計算，并與風洞試驗結果進行對比，如圖3所示。

由圖3可知，升力系數和扭矩系數與風攻角正相關，而阻力系數在風攻角為負時，隨著風攻角的增大而降低，反之增加。受邊界層網格劃分和風特性參數設置等因素影響，模擬所得阻力系數的絕對誤差相對較大，但各工況下三分力系數模擬值與風洞試驗值的相對誤差均不超過20%，且變化趨勢大體相同，表明該模型可進一步用于后續抖振響應數值模擬。

2 抖振響應數值模擬流程與參數設置

2.1 模擬流程

對于二自由度主梁斷面，考慮豎向與扭轉自由度的動力學方程可寫為［10］：

式中 h與α分別表示豎向與扭轉位移；x˙表示變量x對時間求導；mh為單位長度主梁質量；mα為單位長度主梁轉動慣量；ci=2miξiωi （i=h，α）為阻尼系數，其中ωh=1.1567為實測主梁一階正對稱豎彎頻率，ωα=25.244為實測主梁一階正對稱扭轉頻率；ξi（i=h，α）表示阻尼比，各方向均取0.4%；ki=miω2i （i=h，α）為剛度系數；L與M分別表示作用于主梁的氣動升力與扭矩，其線性非定常模型表達如下［10］：

式中 χiu與χiw （i=1，2，3，4）為氣動導納函數，本文偏于保守的取1；K=ωB/U表示折減頻率；H*i與A*i （i=1，2，3，4）表示顫振導數；u表示順風向脈動風速；w表示豎向脈動風速；U表示平均風速；ρ表示空氣密度；B表示梁寬。蘇通大橋主梁斷面抖振響應數值模擬在Fluent中直接完成，流體方程與結構動力學方程耦合求解流程如圖4所示。

圖4中，首先通過NSRFG法生成LES入口湍流，再編制UDF提取每一時間步上作用于主梁斷面的氣動升力與扭矩，并使用4階Runge?Kutta法求解式（1）中主梁斷面動力方程，從而確定主梁斷面的運動狀態。最后，根據Fluent中自帶的Remeshing與Smoothing方法實現網格更新。上述耦合過程由CFD求解器與UDF_BUFFETING間的數據交換體現，屬于雙向弱耦合。

2.2 模擬參數設置

為使數值模型滿足相似準則，設模型幾何縮尺比n=50，則模型中其他參數的相似關系與取值如表2所示［1，16］。LES入口湍流風特性定義如表3所示［14］。

開展抖振數值模擬時，取表3中α=0（常數風剖面），湍流強度則參考蘇通大橋橋址區長期風特性實測結果取值（10 min基本時距湍流強度平均值0.2051）［17］，風場湍流度較高。模擬在配置為96核Intel（R） Xeon（R） Platinum 8163 CPU @ 2.50 GHz 的工作站上完成，采用32核并行計算。首先，關閉動網格，采用LES進行瞬態分析使流場發展穩定。再打開動網格，進行主梁斷面抖振響應FSI模擬，LES時間步長為0.002 s，時間步總數28000。設置變量與連續性方程殘差小于10-5為收斂標準，模擬總耗時超1800 h（2個半月）。LES離散與插值格式設置如表4所示。

3 模擬結果分析與驗證

為開展模擬結果的量化分析與驗證，在主梁節段附近設置監測點BP（-0.5 m， 0 m， 0.16 m），測點處順風向與豎向脈動風速時程如圖5所示，相對應的湍流譜密度如圖6所示。

由圖5和6可知，隨著湍流在流域內不斷發展，模型位置處功率譜密度監測值僅在低頻區與目標值吻合良好，而在高頻區迅速衰減。該現象一方面與數值算法產生的湍動能耗散有關；另一方面由于LES采用濾波函數過濾了風速時程中的高頻成分，從而使LES得到的湍流功率譜在高頻部分偏低。盡管如此，主梁節段模型豎彎與扭轉頻率仍位于LES解析尺度內，故生成的湍流場滿足模型抖振模擬要求。將BP測點脈動風速時程代入式（2），并采用4階Runge?Kutta法求解式（1），獲得考慮與不考慮氣動自激力時主梁節段抖振響應的理論解，其中模型氣動參數根據文獻［18］取值。將其與數值解進行對比，以驗證模擬結果的合理性，如圖7所示。相對應的位移響應功率譜密度分析結果如圖8所示。

由圖7和8可知，當模擬風場湍流度較高時，主梁節段模型風振響應兼具強迫振動（隨機性）與自激振動（周期性）特征，但以強迫振動為主。就時程結果而言，由于抖振力模型存在假設簡化，且數值模型氣動參數與風洞試驗結果存在差異，故數值解與理論解并不完全吻合。當不考慮自激力時（即不考慮流固耦合效應），豎向位移理論解明顯偏大，而扭轉位移理論解則無明顯變化。此外，抖振響應數值解與理論解所包含的頻率成分相近，且主要頻率成分與主梁節段豎彎及扭轉頻率相對應。對隨機振動而言，統計量較時程結果更具代表性，故分別計算了抖振響應數值解與理論解的均方根值（Root Mean Square， RMS），結果如表5所示。

由表5可知，不考慮自激力時的豎向位移RMS相對誤差高達60.7%，而考慮自激力時，豎向位移RMS相對誤差僅1.7%。扭轉位移RMS較小，是否考慮自激力對扭轉位移RMS的影響也較小。實際工程中不考慮自激力獲得的抖振響應總體上偏于安全。

為進一步探究抖振發生時的繞流場細節，以計算域z=160 mm平面（Plane 1）為監測面，任意選取計算穩定時主梁斷面附近瞬時速度場進行展示，并基于Q準則對主梁斷面附近的渦核進行識別，速度場與瞬時渦核分布如圖9所示。其中，位移和加速度的數值模擬結果（縮尺后結果）也一并給出，并將對應時刻的位移結果標注于圖中。

由圖9可知，模型上方渦旋形成于迎風側橋面端點處，而下方渦旋則形成于迎風側風嘴端點處。模型上、下方渦核經充分發展后，均在模型中心線附近達到最大尺寸。特別地，當模型下方渦核破碎消散后，又將在底板與背風側斜腹板交點處重新形成，并向尾渦區發展。此外，模型上、下方渦核發展交替進行。當上方渦核不斷發展時，下方渦核則逐漸破碎消散，反之亦然。而模型附近渦核演化（特征湍流引起）與模型抖振響應間未見明顯關系，表明特征湍流對模型抖振響應的影響較小。

4 結 論

本文基于NSRFG法實現了LES入口湍流生成，通過編制UDF實現了大跨索承橋梁流線型鋼箱梁抖振響應數值模擬，并將抖振響應數值解與理論解進行了對比分析，得到主要結論如下：

（1）三分力系數的CFD計算值與風洞試驗值相對誤差不超過20%，且變化趨勢大體相同，考慮流固耦合效應的抖振響應理論解與數值解的RMS較為接近。

（2）是否考慮自激力對扭轉位移RMS的影響都較小，但不考慮自激力時的豎向抖振位移的RMS則明顯偏大，實際工程中不考慮流固耦合效應獲得的抖振響應總體上偏于安全。

（3）當模擬風場湍流度較高時，主梁節段模型風振響應兼具強迫振動（隨機性）與自激振動（周期性）特征，但以強迫振動為主。模型上、下方渦旋分別形成于迎風側橋面端點和風嘴端點處，兩側渦核經過充分發展后，均在模型中心線附近達到最大尺寸，且下方渦核破碎消散后，又將在底板與背風側斜腹板交點處重新形成，并向尾渦區發展。

（4）模型上、下方渦核發展交替進行。當上方渦核不斷發展時，下方渦核則逐漸破碎消散，反之亦然。模型附近渦核演化（特征湍流引起）與模型抖振響應間未見明顯關系，特征湍流對模型抖振響應的影響較小。

與理論分析及風洞試驗相比，抖振響應數值模擬充分捕捉了模型的繞流場特性，并清晰展示了抖振發生時模型繞流場的演化規律，研究成果有望進一步用于主梁斷面抖振性能優化及橋面風環境評估等。然而，受計算資源限制，數值模擬時的LES解析尺度尚無法完全涵蓋實際風場中的湍流尺度，模擬所得數值風場與實際風場仍存在一定差異。因而，在計算資源允許時，仍可開展更加精細化的數值模擬，以進一步驗證本文結果。

參考文獻

1陳政清. 橋梁風工程［M］. 北京： 人民交通出版社， 2005.

2Ge Y， Xia J， Zhao L， et al. Full aeroelastic model testing for examining wind-induced vibration of a 5，000 m spanned suspension bridge［J］. Frontiers in Built Environment， 2018， 4： 20.

3Simiu E， Scanlan R H. Wind Effects on Structures［M］. New York： John Wiley and Sons， 1996.

4葛耀君. 大跨度橋梁抗風的技術挑戰與精細化研究［J］. 工程力學， 2011， 28（增刊 Ⅱ）： 11-23.

Ge Y J. Technical challenges and refinement research on wind resistance of long-span bridges［J］. Engineering Mechanics， 2011， 28（Sup Ⅱ）： 11-23.

5Ahsan K. Advanced Structural Wind Engineering［M］. New York： Springer， 2013.

6Tamura T. Reliability on CFD estimation for wind-structure interaction problems［J］. Journal of Wind Engineering and Industrial Aerodynamics， 1999， 81（1-3）： 117-143.

7Nariman N A. Influence of fluid-structure interaction on vortex induced vibration and lock-in phenomena in long span bridges［J］. Frontiers of Structural and Civil Engineering， 2016， 10（4）： 363-384.

8詹昊， 廖海黎. 橋梁主梁端部翼板顫振主動控制流固耦合計算［J］. 振動工程學報， 2018， 31（2）： 276-282.

Zhan H， Liao H L. Flutter active control studies of bridge with flaps attached to deck edge by FSI calculation［J］. Journal of Vibration Engineering， 2018， 31（2）： 276-282.

9董國朝， 陳政清， 羅建輝， 等. 安裝亮化燈具導致的斜拉橋拉索風致馳振流固耦合分析［J］. 中國公路學報， 2012， 25（1）： 67-75.

Dong G C， Chen Z Q， Luo J H， et al. Fluid-structure interaction analysis of wind-induced galloping of cables with lamps of cable-stayed bridge［J］. China Journal of Highway and Transport， 2012， 25（1）： 67-75.

10Kavrakov I， Morgenthal G. A synergistic study of a CFD and semi-analytical models for aeroelastic analysis of bridges in turbulent wind conditions［J］. Journal of Fluids and Structures， 2018， 82： 59-85.

11周桐， 楊慶山， 閆渤文， 等. 大氣邊界層大渦模擬入口湍流生成方法綜述［J］. 工程力學， 2020， 37（5）： 15-25.

Zhou T， Yang Q S， Yan B W， et al. Review of inflow turbulence generation methods with large eddy simulation for atmospheric boundary layer［J］. Engineering Mechanics， 2020， 37（5）： 15-25.

12劉震卿， 張沖， 吳曉波， 等. 考慮粗糙地表的近地層風場仿真模擬研究及誤差分析［J］. 空氣動力學學報， 2020， 38（4）： 715-722.

Liu Z Q， Zhang C， Wu X B， et al. Simulation and error analysis of wind field considering surface roughness［J］. ACTA Aerodynamic SINICA， 2020， 38（4）： 715-722.

13Yan B W， Li Q S. Inflow turbulence generation methods with large eddy simulation for wind effects on tall buildings［J］. Computers and Fluids， 2015， 116： 158-175.

14Yu Y， Yang Y， Xie Z. A new inflow turbulence generator for large eddy simulation evaluation of wind effects on a standard high-rise building［J］. Building and Environment， 2018， 138： 300-313.

15Grinderslev C， Lubek M， Zhang Z. Nonlinear fluid-structure interaction of bridge deck： CFD analysis and semi-analytical modeling［J］. Wind and Structures， 2018， 27（6）： 381-397.

16許福友， 馬如進， 陳艾榮， 等. 蘇通大橋全橋氣彈模型設計與模態調試［J］. 工程力學， 2009， 26（12）： 150-154.

Xu F Y， Ma R J， Chen A R， et al. Full aeroelastic model design and modal test for Sutong Bridge［J］. Engineering Mechanics， 2009，26 （12）： 150-154.

17Wang H， Li A， Niu J， et al. Long-term monitoring of wind characteristics at Sutong Bridge site［J］. Journal of Wind Engineering and Industrial Aerodynamics， 2013， 115： 39-47.

18Wang H， Hu R， Xie J， et al. Comparative study on buffeting performance of Sutong Bridge based on design and measured spectrum［J］. Journal of Bridge Engineering， 2013， 18（7）： 587-600.

Buffeting numerical simulation of streamlined steel box girder of long-span cable-supported bridge using fluid-structural interaction

XU Zi-dong 1，2 ?WANG Hao 2 ?LIU Zhen-qing 3

1. Key Laboratory of Concrete and Prestressed Concrete Structures of Ministry of Education， Southeast University， Nanjing 211189，China；

2. School of Civil Engineering， Inner Mongolia University of Science and Technology， Baotou 014010，China；

3. School of Civil and Hydraulic Engineering， Huazhong University of Science and Technology， Wuhan 430074，China

Abstract Evaluation on the buffeting performance of the long-span cable-supported bridge is the significant step for the bridge wind-resistant safety design. However， current buffeting theories are more descriptive， but not explanatory in the sense of revealing basic buffeting physical causes. In addition， the details of the fluid-structure interaction （FSI） during the buffeting are hard to be reproduced. In this work， the Sutong Bridge is selected as the research object， the buffeting responses of its numerical model are simulated in the software Fluent. And the numerical buffeting responses are compared with the theoretical ones. Results show that the forced vibration is the major component of the wind-induced vibration in the wind field with high turbulence intensity. The theoretical and numerical RMS of the buffeting responses are close when considering FSI effect. Vortex core development conducts alternatively between the upper and bottom of the deck. Whats more， signature turbulence has little effect on the buffeting responses of the model.

Keywords long-span bridge; buffeting response; computational fluid dynamics; fluid-structural interaction; large eddy simulation