風速對水線及斜拉索振動的影響

王 劍，畢繼紅，何旭輝，關 健，喬浩玥

(1. 天津城建大學 天津市土木建筑結構防護與加固重點實驗室，天津 300384；2.天津大學 建筑工程學院, 天津 300072；3.中南大學 土木工程學院橋梁工程系，長沙 410000)

隨著現代斜拉橋跨徑的不斷增大，斜拉索越來越長、剛度和阻尼越來越小，導致其在特定風速及降雨條件下產生風雨激振現象。自1988年Hikami等[1]首次觀測到風雨激振現象以來，研究人員對風雨激振的機理進行了大量的研究，并認識到拉索表面水線的形成和振蕩是風雨激振的主要原因。對此Lemaitre等[2-3]等將滑移理論應用到風雨激振研究中，用風壓力和風摩擦力來表示氣流對水膜的作用，模擬水平靜止拉索表面上的水線形成。許林汕等[4]考慮拉索運動對水膜形態變化的影響，根據已有研究成果假設拉索的運動規律，建立液固耦合模型研究水線的運動規律。Taylor等[5]則運用Lemaitre模型，通過數值求解得到隨時間變化的風壓力系數和風摩擦力因數，建立氣液耦合模型分析水平靜止的拉索表面上的水線運動規律。畢繼紅等[6-7]等將滑移理論和單質點運動理論結合，建立了氣液固三相耦合的風雨激振模型。

國內外研究人員通過一系列的現場觀測和風洞試驗研究，發現在降雨強度為小到中雨的條件下，斜拉索僅在一定的風速范圍內發生風雨激振現象，而風速過大或過小均不會出現，即風雨激振有“限速”的特點?，F階段關于風雨激振的起振風速區間的研究大多以風洞試驗為主[8-9]，理論分析方面的研究較少，特別是從水線運動角度的研究就更少。對此，本文對氣液固三相耦合模型的數值計算過程進行改進，采用有限差分法求解水膜運動方程，并且應用有限元軟件COMSOL求解隨時間變化的風壓力系數和風摩擦力因數，在保證計算精度的同時大幅度提高計算效率；在此基礎上，研究風速變化對水線運動及斜拉索振動的影響，分析風雨激振的產生機理。

1 模 型

圖1所示的半徑為R、傾角為α(0°≤α≤90°的斜拉索，受重力和水平方向風的共同作用，風速大小為U0，風偏角為β(0°≤β≤90°)。

圖1 斜拉索空間位置Fig.1 Spatial position of stay cable

取斜拉索任一橫截面A-A為研究對象，如圖2所示。對拉索的二維模型進行研究，忽略重力在拉索軸向方向的影響。

圖2 斜拉索表面水膜受力Fig. 2 Force of water film around cable

根據滑移理論，假設拉索表面存在一層連續的水膜，兩自由度無量綱化的水膜運動方程可表示為

(1)

其中，

(2)

式中：t為時間；μ為水的動力黏度系數；γ為水在空氣中的表面張力系數；h為水膜厚度；Rc為拉索半徑；M為單位長度拉索的質量；ρg和ρ分別為空氣和水的密度；κ為水膜表面的曲率；n和t分別為水膜與空氣交界處的法向向量和切向向量；σg分別為空氣的應力張量；pg為水膜表面所受的空氣壓力。

斜拉索在橫風向和順風向的振動方程為

(3a)

(3b)

其中，

(4a)

(4b)

2 數值求解

2.1 有限差分格式

(5)

其中，

(6a)

(6b)

(6c)

2.2 風壓力系數Cp與風摩擦力因數Cf

水膜運動方程中的風壓力系數Cp和風摩擦力因數Cf是影響水膜運動的關鍵參數。在已有的研究中，畢繼紅等人采用Fluent軟件獲取表面帶有不同形狀水膜的斜拉索的Cp和Cf，計算時間很長，效率極低，僅能獲得前20 s的運動狀態；而Gu等[10]的試驗數據表明，從試驗開始到風雨激振真正形成，需要數十秒甚至上百秒的時間。對此，本文采用有限元軟件COMSOL計算隨時間變化的風壓力系數和風摩擦力因數，大幅提高了計算效率和計算穩定性。

本文應用COMSOL軟件計算時仍采用與原先相同的假設，即假設水膜所受的氣流作用與干燥拉索表面所受的氣流作用相同。每一時間步內均根據上一步計算得到的水膜形狀，將水膜考慮為固態，建立流場中的繞流障礙物，應用Spalart-Allmaras湍流模型，采用穩態計算此時刻的Cp和Cf。

設干燥拉索直徑為d，計算區域及網格劃分如圖3所示。邊界條件設置如下：

左側進口處采用速度邊界條件(u=UN,v=0)；右側出口處采用壓強邊界條件(p=0)；上下邊界處采用完全滑移邊界條件；圓柱表面處采用無滑移邊界條件(u=0,v=0)。

圖3 流場網格劃分Fig.3 Numerical grids of the flow field

將數值求解得到的風壓力系數Cp和風摩擦力因數Cf與Fage等[11]的試驗數據、Celika等[12]及Bi等的數值計算結果進行對比，考慮圓柱直徑為0.16 m,風速為10 m/s的試驗工況，如圖4所示，應用COMSOL軟件求解得出的結果與前人的研究結果非常接近，驗證了數值方法的準確性。

圖4 圓柱表面風壓力系數Cp和風摩擦力因數Cf的分布圖Fig.4 Distribution of wind pressure coefficient Cp and friction coefficient Cf around cylinder

2.3 基本參數及數值計算工況基本

借鑒文獻[13-14]中的風洞試驗數據，基本參數選取如下：斜拉索半徑0.05 m、傾角30°、自振頻率0.952 Hz、線密度8.57 kg/m、阻尼比0.17%；拉索的初始位移、速度和加速度均為零；水膜的初始厚度為0.2 mm；水密度1.0×103kg/m3、水的運動黏性系數1.0×10-6m2/s、水在空氣中的表面張力系數7.2×10-2N/m、空氣密度1.225 kg/m3、空氣的運動黏性系數1.51×10-5m2/s；風偏角22.5°，風速6.0 m/s，6.76 m/s，7.4 m/s，7.5 m/s，7.6 m/s，7.72 m/s，8.2 m/s，8.5 m/s，9.0 m/s，9.5 m/s，10.5 m/s，12.0 m/s，13.5 m/s。

參照已有的研究成果，在采用有限差分法求解水膜運動方程時，時間步長dt=1×10-5s，為保證拉索表面水膜的質量守恒，根據試驗觀測數據，假設水膜的最大厚度為1 mm；為使水膜在拉索表面連續分布，假設水膜的最小厚度為0.02 mm。

2.4 數值求解流程

水膜運動方程(見式(1))為四階偏微分方程，拉索振動方程(見式(3))為二階微分方程，采用差分法求解方程組，應用MATLAB軟件進行數值計算并利用有限元軟件COMSOL求解各個時間步不同水膜形態下的風壓力系數Cp和風摩擦力因數Cf，基本流程如圖5所示。

圖5 數值計算流程圖Fig.5 Flowchart of numerical calculation

3 數值計算結果

3.1 拉索振動響應

圖6顯示了各個風速工況下斜拉索的橫風向及順風向振幅變化。數值計算得到的起振風速區間是6.76 m/s

圖6 風速對拉索振動的影響Fig.6 Influence of wind velocity on cable vibrations

圖7為風速7.5 m/s下拉索在0～150 s內的振動響應。在0～29 s以內，拉索橫風向振幅迅速增大，振幅達到0.19 m左右，如圖6(a)所示；而后振幅迅速衰減，至50 s附近，振幅減至0.13 m；50～100 s內拉索振幅繼續小幅減小，在100 s以后逐步趨于穩定，保持在0.10 m左右，與實驗結果(0.09 m附近)相近，如圖6所示。

7.5 m/s時拉索的順風向振幅如圖7(b)所示。初始時刻，拉索的順風向振幅為0.14 m；前45 s內，拉索振幅迅速減小，最小值為0.03 m；之后小幅增大，72 s后逐漸穩定在0.05 m左右，近似等于拉索橫風向振幅的一半，與現場觀測結果相同。與前人的研究結論相一致，拉索在兩個方向的振動均需要幾十秒(50 s左右)的時間達到穩定振動狀態，因此本文對各風速工況均進行了150 s的計算，選用100～150 s的計算結果進行分析研究。

圖7 風速7.5 m/s時的拉索振動時程曲線Fig.7 Vibration responses of cable with wind velocities 7.5 m/s

3.2 拉索氣動升力和阻力

各個風速工況下的拉索氣動升力及阻力變化幅度曲線，如圖8所示。隨著風速的增加，升力和阻力在總體上均呈現出不斷增大的趨勢。但當風速處于6.76 m/s

對每一風速下拉索的氣動力變化時程進行頻譜分析，如圖9(a)所示。絕大多數風速下升力和阻力的主頻相一致。當風速較低(U0≤6.76 m/s)時，升力和阻力的變化均無明顯的主頻；當風速位于風雨激振的風速區間(6.76 m/s

圖8 風速對氣動力變化最大幅度的影響Fig.8 Influences of wind velocity on aerodynamic forces

圖9(b)為拉索氣動力在其主頻附近的變化幅值與風速變化間的關系曲線。6.76 m/s

綜上所述，決定拉索振動強弱的主要因素是氣動力變化在拉索自振頻率附近幅值，而非其整體的最大變化幅度。

3.3 水線形態及其運動特征

隨著風速從6 m/s增大至9.5 m/s，上水線的最大厚度不斷增加，從0.22 mm逐漸增加至0.81 mm，且其增加的幅度越來越大；至10.5 m/s以上時，上水線的最大厚度達到計算假設的1.0 mm，如圖10所示。上述趨勢與拉索氣動力最大幅度的變化趨勢相近，說明上水線厚度對升力和阻力的變化幅度有顯著影響。

從圖11所示的各風速工況下的上水線振蕩平衡位置可以看出，風速為6 m/s時的上水線中心位置大致在40°附近；隨著風速的增大，上水線平衡位置不斷向背風側偏移；到13.5 m/s風速時移至約70°。文獻[13]和文獻[14]中的坐標系與本文數值計算時選取的坐標系相差90°，經坐標變換后進行對比，發現數值計算結果與風洞試驗結果基本一致，如圖11所示。

圖9 不同風速拉索的氣動力變化頻譜分析Fig.9 Frequency spectrum analysis of aerodynamic forces with various wind velocities

圖10 各風速工況下的上水線厚度Fig. 10 Upper rivulet thickness with various wind velocities

圖11 各風速工況下的上水線振蕩平衡位置Fig.11 Equilibrium positions of upper rivulet oscillations with various wind velocities

圖12所示顯示了各風速工況下的上水線振蕩幅度。由于風速小于7.4 m/s時，上水線振蕩幅度很小，故在此只研究風速大于7.4 m/s時的上水線振蕩幅度。當風速位于風雨激振的風速區間(6.76 m/s

圖12 各風速工況下的上水線振蕩幅度Fig.12 Amplitudes of upper rivulet oscillations with various wind velocities

通過對各風速工況下的上水線位置處水膜厚度變化時程進行頻譜分析(見圖13)，可以進一步研究水線運動與拉索振動間的相互影響。

圖13(a)為上水線振蕩頻率與風速的關系，當風速較低(U0≤6.76 m/s)時，上水線厚度較小且基本無運動；當風速位于風雨激振的風速區間(6.76 m/s

圖13(b)顯示了各風速工況下上水線位置處的水膜厚度變化在其主頻附近的幅值。當風速位于風雨激振的風速區間(6.76 m/s

對比圖6所示的拉索振動與圖10、圖13所示的水線振蕩，可以看出影響拉索振動強弱的主要因素是水線的振動頻率而非其厚度，只有當水線振蕩頻率與拉索自振頻率相近時，才會使得拉索氣動力的頻率接近拉索自振頻率，產生風雨激振現象；此時，拉索的振動幅度與拉索自振頻率附近上水線振蕩的強弱有密切聯系。

分別選取風速6.0 m/s、7.72 m/s、9.0 m/s和13.5 m/s的水膜厚度時程變化用來分析研究水線運動與拉索振動間的相互影響，如圖14所示。

圖13 不同風速拉索的上水線頻譜分析Fig.13 Frequency spectrum analyses of upper rivulets with various wind velocities

6.0 m/s風速下，上、下水線均保持穩定狀態，振蕩幅度很??；7.72 m/s風速下，上、下水線均做周期性的環向振蕩，周期大致與拉索的自振周期一致，且上水線的振蕩幅度遠大于下水線；9.0 m/s風速下，上、下水線亦在振蕩，但其周期顯著減小，遠小于拉索自振周期；至13.5 m/s風速時，上水線振蕩的周期性已不再明顯。

圖15顯示了風速6.0 m/s，7.72 m/s，9.0 m/s和13.5 m/s時上水線位置處水膜厚度變化的頻譜分析。

風速6.0 m/s時，由于上水線基本保持穩定狀態而無明顯的振蕩頻率，導致拉索氣動力變化幅度很小且無明顯主頻，因此拉索振動很弱，不會產生風雨激振現象；風速7.72 m/s時，上水線的環向振蕩顯現出明顯的周期性特征，主頻為0.98 Hz，接近拉索的自振頻率，導致拉索氣動力亦發生同頻率的周期性變化，從而引發拉索的大幅度振動；當風速為9.0 m/s時，上水線的振蕩幅度明顯增大，導致氣動力變化的幅度亦明顯增大，但水線振蕩主頻為2.82 Hz，遠大于拉索自振頻率，致使拉索氣動力變化頻率亦遠大于拉索自振頻率，拉索振幅明顯減??；當風速增大至13.5 m/s時，上水線的振蕩頻率已增大至49.98 Hz，遠遠大于拉索自振頻率，沒有出現風雨激振現象。因此，只有當一定風速范圍內，上水線振蕩頻率與拉索自振頻率相近時，才會導致拉索發生顯著的風雨激振現象。

圖14 拉索表面水膜厚度時程變化Fig.14 Temporal evolution of water film thickness

圖15 上水線頻譜分析Fig. 15 Frequency spectrum analysis of upper rivulet

4 結 論

本文考慮拉索橫風向和順風向的振動對水膜形態變化的影響，將滑移理論和兩自由度振動理論相結合，建立風雨激振時拉索和水膜耦合的理論模型；通過數值計算得到了不同風速工況下150 s的拉索振動響應、氣動力變化和水線運動情況，計算結果與人工降雨條件下拉索節段模型的風洞試驗數據吻合；在此基礎上分析風速對水線運動及拉索振動的影響，得出以下結論：

(1) 由兩自由度模型得到的拉索振動響應、氣動力及水膜形態變化規律均與風洞試驗相近，驗證了模型及數值模擬方法的準確性。

(2) 150 s的拉索振動時程表明拉索在兩個方向的振動均需要經過幾十秒的不穩定期后才會產生風雨激振現象，與已有的現場觀測結論相一致，顯示出提高數值計算效率的必要性。

(3) 隨著風速的增大，水線厚度逐漸增加，上水線環向振蕩的中心位置逐漸向背風側偏移，但拉索只在一定風速范圍內發生大幅度振動，且其橫風向振幅接近順風向振幅的兩倍。

(4) 只有當風速處于特定范圍內時，上水線在拉索表面環向振蕩的幅度最大，且其周期性振蕩的頻率接近拉索自振頻率，導致升力、阻力發生相近頻率的周期性變化，引發拉索在橫風向和順風向發生大幅度振動；而當風速過大或過小時，水線運動及升力、阻力變化的頻率與拉索自振頻率相去甚遠，導致拉索振動明顯減弱。因此，上水線與拉索之間的共振是導致拉索發生風雨激振的主要原因之一。