箱梁顫振過程中旋渦的演化規律及氣動力的變化特性

劉祖軍，楊詠昕，葛耀君

(1.華北水利水電學院土木與交通學院，鄭州 450011；2.同濟大學橋梁工程系，上海 200092)

0 引言

隨著計算方法和計算機性能的不斷提高，計算流體動力學(CFD)技術在橋梁抗風研究中的應用越來越廣泛。Larsen和 Walther采用基于離散渦法模擬了大海帶東橋的風振現象，計算結果與風洞試驗基本吻合［1］。Lucor模擬了低雷諾數下的拉索和主梁的渦振［2］；Lee Sangsan利用CFD方法對韓國的Namhae橋和Seohae橋進行了二維模型分析，所得結果與風洞試驗相近［3］。CFD方法可以獲得十分詳細的流場和壓力分布隨時間變化的信息，研究引起橋梁風振的流場特征，包括壓力、速度分布和旋渦的生成、運動規律，并且在流場顯示圖中可以觀察到旋渦脫落形式，有利于分析流場和結構之間的相互作用。但是CFD數值模擬缺少普遍適用的湍流模型，計算的收斂性和精度還有待提高，因此CFD計算結果的可靠性和準確性需要經過流場測試結果的驗證。

PIV流場測試技術是20世紀80年代發展起來的流動測試手段［4］，具有瞬時性、全場性、無損性和定量性等特點。它既具備了單點測量技術的精度和分辨率，又能獲得平面流場顯示的整體結構和瞬態圖像。

隨著PIV流場測試技術的日臻成熟，已經被廣泛應用到湍流、鈍體繞流、多相流、機翼翼尖渦和火焰等需要獲得復雜流場信息的研究中，因此將PIV技術應用到橋梁風振的流場作用機理研究中，為定量觀測橋梁斷面振動狀態下的旋渦演化規律提供了有利的工具。

同濟大學的張偉采用PIV技術研究了H型斷面風致振動過程中模型表面旋渦的演化過程以及安裝在斷面上的中央穩定板的高度對模型表面旋渦產生的影響［5］。他還根據PIV試驗分析了導流板對箱型橋梁斷面底部和尾部旋渦的影響，得出一些有參考價值的結論［6］，從旋渦演化的角度對斷面的振動進行了闡述，但他沒有分析旋渦的演化對模型表面壓力的影響。

在鈍體空氣動力學中，POD方法常用來簡化和分析從風洞試驗或其它手段獲得的壓力。POD方法在風工程中廣泛應用于低矮建筑、高層建筑以及橋梁的脈動風壓研究。Armitt［7］首先將該方法應用在冷卻塔脈動風壓的研究中。Kareem及Cermark［8-9］將這種方法用于湍流中的矩形截面以及模擬的大氣邊界層中的矩形高層建筑物動風壓分析。同濟大學的周志勇采用POD方法對大跨度橋梁的顫振進行了探討［10］。其研究方法具有啟發意義，但是沒有將壓力場的波動特性與流場變化聯系起來。

扁平箱梁具有較好的流線形，是目前大跨度橋梁主梁設計中采用較多的一種截面形式。主跨1624m的丹麥大海帶東橋是目前跨度最大的單箱截面橋梁，在施工過程中曾出現過超出風洞試驗預期的大幅渦激振動，使得該橋氣動性能的研究倍受關注。

以大海帶東橋的主梁斷面形式為原型，采用PIV流場測試技術并結合CFD數值分析方法，研究了扁平箱梁顫振過程的流場作用機理，應用POD分解技術研究顫振過程中模型表面壓力模態特征函數的分布特性，分析了旋渦演化過程對模型表面氣動力分布特性的影響，以及由此造成的扁平箱梁風致振動的流場作用特征。

1 閉口單箱測振PIV試驗

1.1 PIV試驗設備簡介

PIV試驗在同濟大學土木工程防災國家重點實驗室TJ-4邊界層風洞中進行(圖1，2)，該風洞配套的PIV設備整體引進自美國TSI公司。PIV設備共分為5個部分，即激光器(Dual YAG Laser)、同步儀(Synchronizer)、圖形采集卡及計算機、CCD相機以及粒子投放設備，試驗采用的激光是波長為532nm的綠光。CCD相機采用TSI公司生產的Model 630057 PowerViewTM Plus 2MP PIV相機，最大分辨率為1600×1200像素，每個像素的尺寸是7.4μm×7.4μm，CCD的尺寸為11.8mm×8.9mm。CCD相機獲取圖像對的時間間隔最小為200ns，采樣頻率為15Hz(15pairs/s)。使用的粒子生成液體由丙三醇、丙二醇等原料制成，生成的粒子直徑小于10μm。

圖1 TJ-4風洞的PIV設備示意Fig.1 PIV equipment in the TJ-4 wind tunnel

圖2 TJ-4風洞PIV設備全貌Fig.2 The whole picture of PIV equipment in the TJ-4 wind tunnel

1.2 相位平均方法

在模型振動的任何時刻，其周圍的流場是非定常的，因此對相同振動狀態的流場圖像采用相位平均的處理方法，不但可以消除噪音信號的干擾，而且可以較準確地描述該振動狀態下模型周圍的流場。相位平均方法由Reynolds和Hussain在1972年提出，用來研究流固耦合系統中流體的周期性或準周期性運動［11］。在相位平均方法中，需要一個參考信號用來判斷流體周期性運動的相位，采用激光位移計測得的模型振動位移信號作為參考信號。

1.3 PIV測振試驗及結果分析

箱梁斷面以大海帶東橋的主橋斷面為原型進行縮尺，模型的幾何縮尺比1∶120，頻率比40∶1，試驗的風速比1∶3，風洞流場的湍流度≤0.5%，將模型懸掛在均勻流中進行PIV測振試驗。風洞試驗時風速從零開始直到模型振動發散，試驗過程中雷諾數的變化范圍是0＜Re≤4.8×104，而實橋在實際風場中的雷諾數接近于Re=1.7×106，雖然在一般風洞試驗中，雷諾數相似準則在風洞試驗氣彈相似準則中很重要，但是由于風洞洞體尺寸和試驗風速的限制，風洞試驗雷諾數比實際結構的雷諾數往往要小2～3個數量級。模型截面尺寸如圖3所示，縱向長度0.8m，采用鋁質芯棒，外覆泡沫制成。模型的基本參數如下。

圖3 單箱模型(單位：mm)Fig.3 The model of box(unit：mm)

每延米的質量m=1.3625kg/m，每延米的質量慣性模型矩Im=0.01277kg·m2/m，豎向頻率fh=2.825Hz，扭轉頻率fa=7.451Hz。測振時采用激光位移計記錄模型試驗風速下的位移信號。PIV試驗中相機曝光時間間隔為28μs，激光脈沖延遲時間0.23ms，采樣頻率為15Hz，即每秒獲得15對圖像。在單箱的測振試驗中，記錄了各級風速下的結構位移響應，其結果見圖4。在風速小于18m/s時結構的豎向位移振幅很小，但是存在較小幅度的扭轉振動。當風速超過18m/s后結構振動位移極大值和方差就突然增大，豎向振動的參與程度也不斷增加，直到風速到達20.8m/s時結構出現了顫振失穩現象。

圖5展示了較低風速10m/s時單箱尾端風嘴附近的瞬態流場隨時間的演化過程，每幅圖片的時間間隔為0.067s。從流跡上看尾端風嘴的下側產生尺度較大的橢圓形旋渦，隨后在風嘴的尖角處形成了一個渦量相反的旋渦，結構處于上下側兩個旋渦的共同作用下，接著上側的旋渦消失，下側的旋渦得到了充分的發展，然后又產生了上側旋渦，結構受到上下旋渦的交替控制。但是從試驗結果上看在該風速下風嘴下部旋渦處于主導作用。

圖4 單箱位移響應-風速曲線Fig.4 Displacements vs wind velocity curve of box section

從單箱振動的位移-風速曲線上可以看出在風速18m/s時結構的豎向振動和扭轉振動響應都增加較大。圖6展示了18m/s時單箱尾部風嘴附近流場的瞬態變化過程，從流場變化過程來看，風嘴尖角處的上側旋渦的尺度明顯增加，表明流場內部的能量在該位置處不斷增加，對結構振動起主控作用的上下側旋渦的能量逐漸匹配，輻射到結構的能量也得到加強，該風速下結構的振動受上下旋渦的交替控制作用，模型振幅明顯增大。

圖5 瞬態流場的演化過程(10m/s)Fig.5 The evolution of tail instantaneous flow field(10m/s)

圖6 瞬態流場的演化過程(18m/s)Fig.6 The evolution of tail instantaneous flow field(18m/s)

采用相位平均的方法［10］對單箱尾部風嘴處的旋渦演化規律進行探討。圖7給出了在顫振風速20.8m/s時單箱振動的四分點相位的流場圖，四分點相位的參照變量是扭轉位移，扭轉位移處于峰值時為3π/2相位，最小值時為π/2相位。達到顫振臨界風速時，單箱尾部風嘴的上下側旋渦又重新成為控制結構振動的主導因素，梁底斜邊的位置處旋渦流線的速度梯度加大。

圖7 單箱斷面旋渦驅動結構運動的流場(顫振風速20.8m/s，PIV)Fig.7 Vortex driven structural movement of box section(flutter critical wind speed 20.8m/s，PIV)

2 數值模擬及壓力場的POD分解

2.1 單箱顫振過程的數值模擬

為彌補PIV試驗無法獲得壓力場的不足，通過CFD數值模擬來獲得單箱顫振時表面壓力的變化情況。

數值模擬時采用商業軟件Fluent提供的RANS方法的k-ωSST兩方程模型［12］，計算域的大小參考了同濟大學土木工程防災國家重點試驗室TJ-4風洞中段試驗端設置，計算域沿流線長度為3m(其中上游1m，下游2m)，橫向寬度為0.8m。計算時壁面附近最小網格尺度為0.0004m，計算域采用分塊結構化網格，網格數量為10.6萬。

根據測振試驗得出的顫振臨界狀態下的風速U=20.8m/s、顫振頻率 f=6.054Hz、豎彎幅值 A=1.82cm和扭轉幅值θ=2.4°等參數設定了模型的正弦運動形式，然后根據Fluent的動網格技術求解振動結構周圍的流場和壓力場的變化。計算參數設置為：動量、湍動能和能量耗散均采用兩階迎風格式進行離散，壓力速度耦合采用SIMPLE算法，求解器采用分離式，計算模式選用兩階隱式。邊界條件設定為：速度入口，入口風速20.8m/s，湍流強度與風洞流場相同取0.5% ，壓力出口，計算域的上端和下端設為對稱邊界條件，表面采用無滑移的壁面條件。數值計算的雷諾數與風洞試驗保持一致，即Re=4.8×104，與實橋在實際風場中的雷諾數相差兩個數量級。

圖8給出了通過數值模擬獲得的箱梁尾部旋渦的四分點相位的流場特征圖，與PIV試驗結果對比分析可以發現，本文數值模擬的結果能夠反映箱梁尾部旋渦的演化規律。

圖8 單箱斷面的旋渦驅動結構運動的流場(顫振風速20.8m/s，數值模擬)Fig.8 Vortex driven structural movement of box section(flutter critical wind speed 20.8，CFD)

2.2 壓力場的POD分析

由數值模擬獲得的顫振臨界風速下單箱表面的壓力場就可以采用POD分解方法研究主導壓力場的主要因素。通過對比分析在臨界風速20.8m/s時模型處于振動和靜止兩種條件下表面壓力分布的不同，解釋單箱風致振動的流場機理。

圖9給出了模型振動狀態下表面壓力的POD分解結果。其中第一階模態代表了平均壓力的分布情況，平均壓力分布除兩端部外其它部位比較均勻，并且基本對稱。從第二階模態開始代表了表面壓力的波動部分。由于POD分解的特征值大小反映了各階模態能量的大小，通過分析知第二階模態能量占總波動能量的比例最大(圖9(d))，因此第二階模態是引起結構顫振的主要因素，并且第二階模態展現的壓力分布是不對稱。

圖9 振動狀態時表面壓力POD分解的模態分布(風速20.8m/s)Fig.9 The surface pressure modal distribution of plate in vibration state by POD decomposition(20.8m/s)

圖10 是模型靜止狀態下表面壓力POD分解后的模態分布情況。第一階模態的分布情況與振動狀態下基本相同，而第二階模態的分布情況在顫振狀態和靜止時有較大的差別(圖9(b)，圖10(b))，在振動狀態時主控波動壓力向迎風端漂移，對這一現象從流場機制可以這樣解釋：氣流流過單箱斷面時在其尾部風嘴附近產生明顯的旋渦，當模型固定時其尾部靜止不動，造成尾端與旋渦的相對運動較大，而模型運動時尾部旋渦基本控制了模型尾部的運動趨勢，因此二者之間的相對運動較小，這樣就減輕了該區域氣體分子無規則運動的程度，降低了該區域流動的紊亂程度，從而使尾部表面的波動壓力減小。

圖10 固定狀態時表面壓力POD分解的模態分布(風速20.8m/s)Fig.10 The surface pressure modal distribution of plate in fixed state by POD decomposition(20.8m/s)

圖11 振動過程中旋渦的推動作用Fig.11 Vortex driven structural moment in the process of model vibration

3 結論

通過對箱梁顫振過程的PIV測振試驗并結合數值模擬計算，采用POD分解技術分析顫振臨界狀態下的模型表面壓力場的分布特性，得出以下幾點結論：

(1)箱梁顫振過程中，模型的振動與尾部旋渦演化存在一定的聯系。當風速較低時，模型的振幅較小，尾部風嘴附近上側旋渦的尺度也很小，而下側旋渦尺度較大，形狀接近于圓形，隨著風速的增加，模型尾部風嘴上側旋渦尺寸明顯增大，達到與下側旋渦尺寸相當的地步，并且上下側的旋渦交替出現，此時結構的振幅顯著增大，上下側旋渦的交替作用主導了結構振動直到振動發散。

(2)箱梁處于顫振臨界狀態時，模型振幅很大，其尾部旋渦的演化具有一定的規律性。在一個完整的振動周期內，尾部旋渦的演化過程為：當處于振動平衡位置時，模型尾部風嘴附近的下側存在著一個接近于圓形的旋渦，而上側的旋渦還沒有出現；在由平衡位置運動到負扭轉最大位移時，模型尾部風嘴下側的旋渦尺度急劇增加，對結構產生較強的向上推動作用，而上側的旋渦在逐漸形成；隨后結構由負最大位移回到平衡位置時，上側旋渦逐漸消失，下側旋渦的尺度也不斷減??；然后模型由平衡位置運動到正的最大位移時，其尾部風嘴附近上下側旋渦的尺寸都不斷增加，并最終達到兩者尺度達相當的程度(圖11)。

(3)壓力場的POD分解表明在箱梁顫振過程中存在著主控波動壓力向迎風側風嘴漂移的特點。當模型靜止時模型尾端的波動壓力較大，而迎風端較小，但是當模型處于顫振臨界狀態時，迎風端的波動壓力值最大，尾部較小，形成主控波動壓力向迎風端的漂移。

(4)箱梁顫振過程中模型表面主控波動壓力向迎風側風嘴漂移的特點也正說明了箱梁迎風側風嘴的形狀對其振動穩定性的影響很大。

［1］ LARSEN A，WALTHER J.Aeroelastic analysis of bridge grider sections based on discrete vortex simulations［J］.Journal of Wind Engineering and Industrial Aerodynamics，1997，67＆68：253-265.

［2］ LUCOR D，IMAX L，KAMIADAKIS G E.Vortex disloca-tions and force distribution of long flexible cylinders subjected to sheared flows［J］.Journal of Fluids and Structures，2001，15(3-4)：641-650.

［3］ LEE SS，LEE JS，KIN J D.Prediction of vortex-induced wind loading on long-span bridges［J］.Journal of Wind Engineering and Industrial Aerodynamics，1997，67＆68：267-278.

［4］ ADRIAN R J.Twenty years of particle image velocimetry［J］.Experiments in Fluids，2005，39(2)：159-169.

［5］ 張偉，葛耀君.H形斷面橋梁風致振動的流場機理［J］.土木工程學報，2009，42(5)：90-95.

［6］ 張偉，葛耀君.導流板對大跨度風振影響的流場機理［J］.中國公路學報，2009，22(3)：52-56.

［7］ ARNLITT J.Eigenvector analysis of pressure fluctuations on the west Burton instrumented cooling tower［R］.Central Electricity Research Laboratories(UK)Internal Report，RD/L/N114/68，1968.

［8］ KAREEM A，CERMAK J E.Pressure fluctuations on a square building model in boundary-layer flows［J］.J.Wind Eng.Ind.Aerodyn.，1984，16：17-41.

［9］ BERKOOZ G，HOLMES P，LUMLEY J L.Proper orthogonal，decomposition in the analysis of turbulent flows［J］.Annual Review of Fluid Mechanics，1993，25：537.

［10］周志勇，鄧世為，楊立坤.基于POD技術對大跨度橋梁顫振機理的探討［C］//第十四屆全國結構風工程學術會議論文集，2009：221-227.

［11］FUJISAWA N，TAKEDA G，IKE N.Phase-averaged characteristics of flow around a circular cylinder under acoustic excitation control［J］.Journal of Fluids and Structures，2004，19：159-170.

［12］ BOURIS D，BERGELES G.2D LES of vortex shedding from a square cylinder［J］.Journal of Wind Engineering and Industrial Aerodynamics，1999，80(1-2)：31-46.