風(fēng)速對水線及斜拉索振動的影響

王 劍，畢繼紅，何旭輝，關(guān) 健，喬浩玥

(1. 天津城建大學(xué) 天津市土木建筑結(jié)構(gòu)防護(hù)與加固重點(diǎn)實(shí)驗(yàn)室，天津 300384；2.天津大學(xué) 建筑工程學(xué)院, 天津 300072；3.中南大學(xué) 土木工程學(xué)院橋梁工程系，長沙 410000)

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

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

1 模 型

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

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

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

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

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

(1)

其中，

(2)

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

斜拉索在橫風(fēng)向和順風(fēng)向的振動方程為

(3a)

(3b)

其中，

(4a)

(4b)

2 數(shù)值求解

2.1 有限差分格式

(5)

其中，

(6a)

(6b)

(6c)

2.2 風(fēng)壓力系數(shù)Cp與風(fēng)摩擦力因數(shù)Cf

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

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

設(shè)干燥拉索直徑為d，計(jì)算區(qū)域及網(wǎng)格劃分如圖3所示。邊界條件設(shè)置如下：

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

圖3 流場網(wǎng)格劃分Fig.3 Numerical grids of the flow field

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

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

2.3 基本參數(shù)及數(shù)值計(jì)算工況基本

借鑒文獻(xiàn)[13-14]中的風(fēng)洞試驗(yàn)數(shù)據(jù)，基本參數(shù)選取如下：斜拉索半徑0.05 m、傾角30°、自振頻率0.952 Hz、線密度8.57 kg/m、阻尼比0.17%；拉索的初始位移、速度和加速度均為零；水膜的初始厚度為0.2 mm；水密度1.0×103kg/m3、水的運(yùn)動黏性系數(shù)1.0×10-6m2/s、水在空氣中的表面張力系數(shù)7.2×10-2N/m、空氣密度1.225 kg/m3、空氣的運(yùn)動黏性系數(shù)1.51×10-5m2/s；風(fēng)偏角22.5°，風(fēng)速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。

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

2.4 數(shù)值求解流程

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

圖5 數(shù)值計(jì)算流程圖Fig.5 Flowchart of numerical calculation

3 數(shù)值計(jì)算結(jié)果

3.1 拉索振動響應(yīng)

圖6顯示了各個(gè)風(fēng)速工況下斜拉索的橫風(fēng)向及順風(fēng)向振幅變化。數(shù)值計(jì)算得到的起振風(fēng)速區(qū)間是6.76 m/s

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

圖7為風(fēng)速7.5 m/s下拉索在0～150 s內(nèi)的振動響應(yīng)。在0～29 s以內(nèi)，拉索橫風(fēng)向振幅迅速增大，振幅達(dá)到0.19 m左右，如圖6(a)所示；而后振幅迅速衰減，至50 s附近，振幅減至0.13 m；50～100 s內(nèi)拉索振幅繼續(xù)小幅減小，在100 s以后逐步趨于穩(wěn)定，保持在0.10 m左右，與實(shí)驗(yàn)結(jié)果(0.09 m附近)相近，如圖6所示。

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

圖7 風(fēng)速7.5 m/s時(shí)的拉索振動時(shí)程曲線Fig.7 Vibration responses of cable with wind velocities 7.5 m/s

3.2 拉索氣動升力和阻力

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

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

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

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

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

3.3 水線形態(tài)及其運(yùn)動特征

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

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

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

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

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

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

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

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

圖13(a)為上水線振蕩頻率與風(fēng)速的關(guān)系，當(dāng)風(fēng)速較低(U0≤6.76 m/s)時(shí)，上水線厚度較小且基本無運(yùn)動；當(dāng)風(fēng)速位于風(fēng)雨激振的風(fēng)速區(qū)間(6.76 m/s

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

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

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

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

6.0 m/s風(fēng)速下，上、下水線均保持穩(wěn)定狀態(tài)，振蕩幅度很小；7.72 m/s風(fēng)速下，上、下水線均做周期性的環(huán)向振蕩，周期大致與拉索的自振周期一致，且上水線的振蕩幅度遠(yuǎn)大于下水線；9.0 m/s風(fēng)速下，上、下水線亦在振蕩，但其周期顯著減小，遠(yuǎn)小于拉索自振周期；至13.5 m/s風(fēng)速時(shí)，上水線振蕩的周期性已不再明顯。

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

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

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

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

4 結(jié) 論

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

(1) 由兩自由度模型得到的拉索振動響應(yīng)、氣動力及水膜形態(tài)變化規(guī)律均與風(fēng)洞試驗(yàn)相近，驗(yàn)證了模型及數(shù)值模擬方法的準(zhǔn)確性。

(2) 150 s的拉索振動時(shí)程表明拉索在兩個(gè)方向的振動均需要經(jīng)過幾十秒的不穩(wěn)定期后才會產(chǎn)生風(fēng)雨激振現(xiàn)象，與已有的現(xiàn)場觀測結(jié)論相一致，顯示出提高數(shù)值計(jì)算效率的必要性。

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

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