王劍 畢繼紅 關(guān)健 喬浩玥 邵倩 周燕 管青海

摘要：基于滑移理論推導(dǎo)三維水膜運動方程，重點考慮斜拉索表面水膜形態(tài)變化對風(fēng)壓力系數(shù)和風(fēng)摩擦力系數(shù)的影響，建立模擬風(fēng)雨條件下三維節(jié)段斜拉索表面水膜形態(tài)變化的理論模型。在此基礎(chǔ)上推導(dǎo)出三維節(jié)段斜拉索氣動力的計算公式，研究水線運動對拉索氣動力及其振動響應(yīng)的影響。結(jié)果表明：應(yīng)用三維滑移理論模型數(shù)值模擬得到的拉索表面上下水線的位置、形態(tài)與風(fēng)洞試驗結(jié)果相近;發(fā)生風(fēng)雨激振時斜拉索的順風(fēng)向振幅明顯小于橫風(fēng)向振幅;水線沿拉索環(huán)向與軸向運動使得水膜形態(tài)發(fā)生低頻周期性變化，導(dǎo)致拉索氣動升力與阻力的同頻周期性變化，從而引起拉索的大幅度振動，產(chǎn)生風(fēng)雨激振現(xiàn)象。

關(guān)鍵詞：風(fēng)雨激振;斜拉索;三維模型;滑移理論;水線

中圖分類號：TU312+。1;U443.38文獻(xiàn)標(biāo)志碼：A 文章編號：1004-4523（2020）03-0559-11

DOI：10.16385/i.cnki.issn.1004-4523.2020.03.015

引言

斜拉橋拉索由于具有質(zhì)量輕、剛度小和阻尼小的特點，在風(fēng)雨共同作用下，極易發(fā)生大幅度的低頻振動，即風(fēng)雨激振現(xiàn)象。國內(nèi)外諸多學(xué)者通過一系列的現(xiàn)場觀測和風(fēng)洞試驗研究，發(fā)現(xiàn)斜拉索表面上水線的形成與振蕩是風(fēng)雨激振現(xiàn)象的主要標(biāo)志和關(guān)鍵影響因素。為此，科研人員一方面采用人工水線風(fēng)洞試驗測量水線對拉索氣動性能的影響，另一方面通過人工降雨風(fēng)洞試驗并采用各種測量方式觀測水線的運動特性。

由于水線尺寸較小，形態(tài)變化多種多樣，僅僅依靠風(fēng)洞試驗研究很難全面掌握水線運動與斜拉索振動間的內(nèi)在聯(lián)系，因此有必要采用數(shù)值模擬方法進(jìn)行相應(yīng)的研究。自2007年Lemaitre等首次應(yīng)用滑移理論研究水平靜止拉索表面上的水膜形態(tài)變化并模擬水線的形成以來，經(jīng)過Taylor等、許林汕等、畢繼紅和王劍等的不斷完善，滑移理論逐漸發(fā)展成為數(shù)值模擬水線運動的重要方法。

目前，滑移理論的研究還主要集中于二維斷面模型研究，而斜拉索的風(fēng)雨激振是復(fù)雜的三維問題，應(yīng)考慮水的軸向流動和氣流的軸向流動。對此，Bi等基于滑移理論推導(dǎo)出了考慮拉索振動和氣流作用的三維節(jié)段拉索表面的水膜運動方程，將拉索振動響應(yīng)作為已知條件帶入水膜運動方程，重點研究拉索振動對水線運動的影響規(guī)律。眾所周知，水線運動受重力、氣流作用和拉索振動三個方面共同作用，而水膜形態(tài)變化又會對氣流作用產(chǎn)生影響。對此，本文應(yīng)用已建立的三維水膜運動方程，忽略拉索振動對水膜形態(tài)的影響，而重點考慮斜拉索表面水膜形態(tài)變化對風(fēng)壓力系數(shù)和風(fēng)摩擦力系數(shù)的影響，建立模擬風(fēng)雨條件下三維節(jié)段斜拉索表面水膜形態(tài)變化的理論模型;在此基礎(chǔ)上推導(dǎo)出三維節(jié)段斜拉索氣動力的計算公式，研究水線運動對拉索升力、阻力的影響;而后將拉索氣動力施加到斜拉索模型上得到斜拉索的振動響應(yīng);與已有的風(fēng)洞試驗及數(shù)值模擬結(jié)果進(jìn)行比對驗證，揭示風(fēng)雨激振現(xiàn)象的產(chǎn)生機理。相比于以往的基于滑移理論的研究，本文突破了二維斷面模型的限制，建立的三維模型能夠考慮水線沿斜拉索的軸向流動，更加接近實際，可靠性更高。

1 模型

1.1水膜運動方程

參考文獻(xiàn)推導(dǎo)水膜運動方程，半徑為R、傾角為a（0°≤a≤90°）的斜拉索節(jié)段，受水平方向氣流和重力的共同作用，如圖1所示，風(fēng)速為U，風(fēng)偏角為β（0°≤β≤90°）。

將重力分解，則作用在斜拉索斷面內(nèi)的重力分量g_N和沿斜拉索軸向的重力分量g_z分別為

采用圖2所示的柱坐標(biāo)系（e_r，e_θ，e_z），根據(jù)滑移理論，假設(shè)斜拉索表面存在一層連續(xù)的水膜，斜拉索節(jié)段任一斷面（A-A）上的水膜受力如圖3所示。水膜內(nèi)任一點的坐標(biāo)為（r，θ，z），R≤r≤R+h，h為水膜厚度;速度u表示為分量形式u=u_rr_r+u_θe_θ十u_ze_r，則三維Navier-Stokcs公式可寫為

方程的邊界條件包括：

（1）位移邊界條件：水膜在底面（r-R）處相對于斜拉索表面靜止，即

式中 t為時間，u為水的動力黏度系數(shù)，r為水在空氣中的表面張力系數(shù)（r，u與θ無關(guān)），δ和δ_g分別為水膜和空氣的應(yīng)力張量，I為單位向量，n為水膜與空氣交界處的法向向量，K為水膜表面的曲率（K=-▽·n），p為水膜內(nèi)的壓強;p_g為水膜表面所受空氣壓力，τ_g為空氣黏滯力張量。

（3）自由邊界條件：

水膜在與空氣交界面F（r，θ，Z，t）=R+h（θ，z，t）-r=0處滿足

將邊界條件帶人式（4）并進(jìn)行無量綱處理，得到無量綱的三維水膜運動方程

1.2 拉索氣動力與拉索振動方程

斜拉索任一斷面（A-A）的受力如圖4所示，圖中F_r（θ，z）和F_θ（θ，Z）分別為水膜底面（r=R）處的法向力和切向力。水膜底面（r=R）處的應(yīng)力張量為：

整個三維剛性斜拉索節(jié)段的氣動阻力和升力分別為

2 數(shù)值求解

2.1 風(fēng)壓力系數(shù)C_p與風(fēng)摩擦力系數(shù)C_f

在滑移理論中，氣流對水膜形態(tài)變化的影響主要體現(xiàn)在水膜運動方程中的風(fēng)壓力系數(shù)C_p和風(fēng)摩擦力系數(shù)C_f;而這兩個參數(shù)又會隨著水膜形態(tài)的變化而改變，如圖5所示。Bi等在應(yīng)用三維模型研究風(fēng)雨激振時忽略了水膜形態(tài)變化對風(fēng)壓力系數(shù)C_p和風(fēng)摩擦力系數(shù)C_f的影響，采用固定參數(shù)，具有一定的局限性，無法體現(xiàn)氣流與水膜運動之間的相互作用。由于風(fēng)雨激振中的水線位置和形狀多種多樣，很難通過風(fēng)洞試驗來確定每一時刻的風(fēng)壓力系數(shù)和風(fēng)摩擦力系數(shù)。因此，本文仍采用二維斷面模型的相關(guān)研究方法，暫時忽略軸向氣流的作用，應(yīng)用COMSOL軟件計算每一斷面上隨時間變化的風(fēng)壓力系數(shù)和風(fēng)摩擦力系數(shù)。具體參數(shù)設(shè)置詳見文獻(xiàn)。

2.2 數(shù)值計算流程

圖6顯示了采用MATLAB軟件數(shù)值求解三維水膜運動方程（式（8））及計算拉索氣動力（式（11））和拉索振動響應(yīng)（式（12））的基本流程。

計算C_p和C_f時在COMSOL軟件中劃分流場網(wǎng)格，此過程需要用到上一時刻的水膜形態(tài)，因此在劃分網(wǎng)格時需進(jìn)行參數(shù)化處理;水膜運動方程（式（8））為四階非線性偏微分方程，采用預(yù)測-校正差分格式數(shù)值求解;依據(jù)式（11a）和（11b），在環(huán)向和軸向采用二重數(shù)值積分方法分別計算拉索的氣動阻力F_x和升力F_v;根據(jù)拉索振動方程（式（12）），采用四階Runge-Kutta法求解順風(fēng)向與橫風(fēng)向振動響應(yīng)。

2.3 基本參數(shù)

Li等采用超聲波測厚系統(tǒng)在人工降雨風(fēng)洞試驗中測量了風(fēng)雨激振時斜拉索表面的水膜厚度分布，為此后數(shù)值研究水膜形態(tài)變化提供了實驗基礎(chǔ)。參照該試驗，本文研究選取的基本參數(shù)如下：斜拉索半徑R=0.05m，自振頻率f₀=0.952Hz，風(fēng)偏角β=22.5°，水膜初始厚度h₀=0.2mm，重力加速度g=9.8m/s²，水密度p=1.0×10³kg/m³，水的運動黏性系數(shù)v=1.0×10^-6m²/s，水在空氣中的表面張力系數(shù)γ=7.2×10^-2N/m，空氣密度p_g=1.225kg/m³，空氣的運動黏性系數(shù)v_g=1.51×10^-5m²/s，斜拉索的線密度p_s=8.57kg/m，阻尼比ξ₀=0.17%。

L1等在試驗中發(fā)現(xiàn)，傾角α=30°、風(fēng)偏角β=22.5‘的斜拉索在風(fēng)速U=6.76-8.04m/s時會發(fā)生風(fēng)雨激振現(xiàn)象，在風(fēng)速U=7.72m/s時拉索振動尤為明顯。故本文在數(shù)值計算過程中選擇傾角（α=30°、風(fēng)偏角β=22.5°、風(fēng)速U=7.72m/s作為計算參數(shù)。

根據(jù)文獻(xiàn)中三維水膜方程的求解精度研究，本文選擇0.5m長的三維剛性節(jié)段斜拉索作為研究對象;沿斜拉索軸向劃分為25個斷面，每個斷面環(huán)向離散為120個點，時間步長為10^-3s;水膜的初始厚度0.1mm。忽略雨滴的影響，參考文獻(xiàn)假設(shè)三維節(jié)段拉索表面的水量守恒，總水量（體積y）不隨時間變化，即滿足

應(yīng)用差分法求解水膜運動方程需要水膜在斜拉索表面連續(xù)分布，為此設(shè)定水膜的最小厚度h_min=0.01mm。

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

以往的研究表明，風(fēng)雨激振的形成需要數(shù)10s以上的時間，故本文進(jìn)行了150s的數(shù)值計算，并選取120-140s內(nèi)的計算結(jié)果進(jìn)行研究。

3.1水線的形成

圖7顯示t=0.9s時拉索表面的水膜厚度分布情況。為便于觀察，將三維斜拉索圓柱體展開，橫軸表示斜拉索軸向方向，縱軸表示環(huán)向方向，圖中顏色表示水膜厚度。t=0.9s時，水膜在拉索下方聚集形成下水線，范圍在274.5°-280°之間，寬度約為8.8mm，厚度約為0.86mm，此時還未形成明顯的上水線。至t=1.4s時，在74°-94°范圍內(nèi)出現(xiàn)了明顯的上水線，寬度約為13.9mm，厚度約為0.33mm，如圖8所示。

3.2 三維節(jié)段斜拉索表面的水膜形態(tài)變化

圖9-14顯示了t=136.0-137.0s時間段內(nèi)斜拉索表面的水膜厚度分布。下水線位于275°-280°之間，最大厚度約為1.48mm。其位置、形態(tài)變化很小，只有厚度在斜拉索軸向方向有明顯變化，體現(xiàn)出了氣動力和重力共同作用下水的軸向流動。上水線則在55°-90°之間運動，最大厚度約為0.46mm，與風(fēng)洞試驗的觀測結(jié)果相近。與試驗結(jié)果不同的是，數(shù)值計算結(jié)果中上水線的環(huán)向運動幅值小于試驗結(jié)果，其原因主要是計算中未考慮拉索振動對水膜形態(tài)變化的影響;此外，數(shù)值計算結(jié)果中還出現(xiàn)了水滴自上水線位置沿迎風(fēng)側(cè)下滑并匯合到下水線的現(xiàn)象。產(chǎn)生這一現(xiàn)象的因素除計算中未考慮拉索振動對水膜形態(tài)變化的影響外;另一因素可能是應(yīng)用滑移理論數(shù)值模擬風(fēng)雨激振時無法模擬降雨過程中雨滴下落至斜拉索表面的現(xiàn)象。

3.3 水線的環(huán)向運動

圖15一18分別表示t=120-140s內(nèi)Z=0.1m，Z=0.2m，Z=0.3m和Z=0.4m斷面處的水膜厚度時程變化。圖像橫軸表示時間，縱軸表示拉索環(huán)向位置，顏色表示水膜厚度?？梢钥吹竭@四幅圖像區(qū)別不大。首先，與試驗結(jié)果相一致，下水線的形態(tài)、厚度幾乎一致，均穩(wěn)定在相同范圍內(nèi)（273°-281°）;其次，上水線寬度明顯大于下水線，而厚度明顯小于下水線，與試驗觀測結(jié)論一致;再次，上水線的運動范圍在40°-90°之間，中心位置與試驗結(jié)果一致而范圍略大;最后，這4個斷面處均有水從上水線沿迎風(fēng)側(cè)下滑匯聚至下水線，與二維斷面研究的結(jié)果一致。

為研究上水線位置處水膜形態(tài)變化的運動特性，分別研究Z=0.1m，Z=0.2m，Z=0.3m和Z=0.4m四個斷面θ=63°位置處的水膜厚度，如圖19-22所示。這四個斷面上水線位置處的水膜厚度變化范圍均在0.1-0.4mm范圍內(nèi);對其分別進(jìn)行頻譜分析，發(fā)現(xiàn)水膜厚度變化的主頻主要有三個：0.253，0.558和0.802Hz，均小于拉索自振頻率（0.952Hz），與試驗結(jié)果有一定區(qū)別。對比文獻(xiàn)僅考慮拉索振動影響的研究，可以發(fā)現(xiàn)產(chǎn)生這一差別的主要原因是本文的計算中未考慮拉索振動的影響，說明拉索振動是導(dǎo)致水線周期性環(huán)向振蕩的主要因素之一。

由于在計算中出現(xiàn)了水從迎風(fēng)側(cè)自上水線下滑至下水線的現(xiàn)象，故在上述四個斷面上的迎風(fēng)側(cè)各選取同一位置（θ=9°）的水膜厚度作為研究對象，如圖23-26所示。這4個位置處的水膜厚度變化范圍在0.1-0.5mm內(nèi)，略大于上水線位置處的水膜厚度變化幅度。頻譜分析顯示迎風(fēng)側(cè)的水膜厚度變化主頻與上水線位置處一致，仍為0.253，0.558和0.802Hz。比較這4個斷面8個位置的水膜厚度變化，可以發(fā)現(xiàn)均有0.802Hz這個頻率。

3.4 水線的軸向運動

為研究水線的軸向運動，分別選取t=127-132s內(nèi)Z=0.30m和Z=0.32m兩個相鄰斷面作為研究對象。由于上水線在40°-90°之間周期性振蕩，同時有水周期性的沿迎風(fēng)側(cè)從上水線下滑匯聚到下水線，故研究上水線（θ=60°和θ=63°）及迎風(fēng)側(cè)（θ=9°和（θ=12°）的水膜厚度時程變化，如圖27所示。Z=0.30m和x=0.32m是數(shù)值計算中的兩個相鄰斷面，圖27中的對比顯示這兩個斷面的相同環(huán)向位置處的水膜厚度變化的幅值基本一致，且迎風(fēng)側(cè)的變化幅度明顯大于上水線位置處。

在同一角度位置，兩條水膜厚度時程曲線相近，但水膜最大厚度出現(xiàn)的時間明顯不同，有明顯的時間差存在，說明水在沿著斜拉索軸向流動;比較同一斷面內(nèi)相鄰角度的時程曲線，可以發(fā)現(xiàn)水膜厚度變化亦非完全同步。因此，通過三維模型的數(shù)值計算可以看出，水線在拉索表面同時沿環(huán)向和軸向運動。事實上，水線沿拉索環(huán)向和軸向運動是在重力、氣流和拉索振動共同作用下形成的，但在本文的研究中，忽略了拉索振動的影響，僅研究重力和氣流作用下的水線運動。

3.5 拉索氣動力

圖28和29分別為t=120-140s內(nèi)拉索氣動升力與阻力的時程變化曲線和頻譜分析。兩者的變化趨勢一致，但升力的變化幅值較大，升力的變化范圍是—0.65-2.45N，阻力的變化范圍是4.50-6.90n.頻譜分析顯示升力與阻力的主頻均為0.253，0.558和0.802Hz，與上水線位置處的水膜厚度變化頻率一致，亦與水線下滑導(dǎo)致的迎風(fēng)側(cè)水膜厚度變化頻率一致。這一現(xiàn)象說明水膜形態(tài)變化是導(dǎo)致拉索氣動力變化的直接重要因素。

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

圖30所示的t=120-140s內(nèi)拉索橫風(fēng)向振動響應(yīng)顯示橫風(fēng)向振幅約為0.114m，與風(fēng)洞試驗觀測數(shù)據(jù)基本一致。對比圖31顯示的順風(fēng)向振動的時程曲線，可以發(fā)現(xiàn)拉索發(fā)生風(fēng)雨激振時的順風(fēng)向振幅約為0.054m，明顯小于橫風(fēng)向振幅（0.114m），與Ni等的現(xiàn)場觀測結(jié)果及以往的二維斷面模型研究結(jié)論相一致，產(chǎn)生這一現(xiàn)象的主要原因是升力的變化幅值明顯大于阻力。由于氣動升力有正有負(fù)，故拉索橫風(fēng)向振動的平衡位置大致位于原點附近的0.012m處;而阻力一直為正，導(dǎo)致順風(fēng)向振動的平衡位置遠(yuǎn)離原點，大致位于0.142m處。

結(jié)合水膜形態(tài)研究，可以看出風(fēng)雨激振是由于水膜形態(tài)（厚度）的低頻周期性變化導(dǎo)致拉索氣動升力與阻力的同頻周期性變化，進(jìn)而引起拉索的大幅度振動。

4 結(jié)論

本文基于滑移理論，忽略拉索振動對水膜形態(tài)的影響，應(yīng)用三維水膜運動方程研究重力和氣流共同作用下三維斜拉索表面的水膜形態(tài)變化;并推導(dǎo)出三維節(jié)段斜拉索氣動力的計算公式，研究水線運動對拉索升力、阻力及拉索振動的影響，得到下列結(jié)論：

（1）應(yīng)用三維滑移理論模型數(shù)值模擬得到的拉索表面上下水線的位置、形態(tài)與實驗結(jié)果相近。

（2）在三維節(jié)段模型中考慮水膜形態(tài)變化對風(fēng)壓力系數(shù)和風(fēng)摩擦力系數(shù)影響，可以較好地模擬水線沿拉索環(huán)向與軸向運動。上水線的環(huán)向振蕩及水線下滑導(dǎo)致上水線附近及迎風(fēng)側(cè)水膜厚度的周期性變化，二者主頻一致，均小于拉索自振頻率。

（3）阻力的變化幅度小于升力，導(dǎo)致斜拉索發(fā)生風(fēng)雨激振時的順風(fēng)向振幅明顯小于橫風(fēng)向振幅，與已有的研究結(jié)果相一致。

（4）拉索氣動升力與阻力的變化頻率與上水線及迎風(fēng)側(cè)位置處水膜厚度變化頻率一致。水膜形態(tài)變化是導(dǎo)致拉索氣動力變化的直接重要因素。水膜形態(tài)周期性變化導(dǎo)致拉索氣動力的周期性變化，引發(fā)拉索橫風(fēng)向與順風(fēng)向的大幅振動，形成風(fēng)雨激振現(xiàn)象。