邵亞會,趙德勝,趙 彤,翟奧博,汪亦顯

(合肥工業(yè)大學(xué) 土木與水利工程學(xué)院,安徽 合肥 230009)

π型梁結(jié)構(gòu)廣泛應(yīng)用于懸索橋建造中,其截面為開口形式,風(fēng)荷載作用下,尾流易產(chǎn)生周期性的渦旋脫落,加之其抗扭剛度差、質(zhì)量輕,易產(chǎn)生較大幅度的渦激振動。目前對π型主梁渦振的相關(guān)研究成果較豐富。文獻(xiàn)[1]整合風(fēng)洞試驗、數(shù)值模擬的數(shù)據(jù)樣本,建立多種神經(jīng)網(wǎng)絡(luò)模型,對π型梁渦激振動性能進(jìn)行智能識別;文獻(xiàn)[2]基于結(jié)構(gòu)設(shè)計,通過調(diào)整π型梁底板之間距離,探究π型梁與空氣動力相關(guān)的特性;文獻(xiàn)[3]以某斜拉橋為例,采用階段模型方法,研究穩(wěn)定板和擾流板對π型截面橋的影響;文獻(xiàn)[4]通過風(fēng)洞試驗發(fā)現(xiàn),扭轉(zhuǎn)-豎彎頻率比接近1時,π型主梁在垂直和扭轉(zhuǎn)方向上表現(xiàn)出耦合振動。

由于現(xiàn)場實測和風(fēng)洞試驗的代價昂貴,數(shù)值模擬是現(xiàn)階段研究大跨度橋梁渦激振動的主要方法之一[5]。數(shù)值模擬用于橋梁渦振模擬及渦振機理研究時,結(jié)果準(zhǔn)確度較高,且具有可視化的優(yōu)點。隨著橋梁跨度的不斷增大,在概率背景下評估風(fēng)振成為設(shè)計特大跨度橋梁需要考慮的問題。文獻(xiàn)[5]結(jié)合Scanlan渦旋脫落半經(jīng)驗?zāi)Ｐ蛯皹虻臏u振發(fā)生概率進(jìn)行評估;文獻(xiàn)[6]使用多種一次二階矩法進(jìn)行顫振失穩(wěn)的可靠性分析;文獻(xiàn)[7]考慮到風(fēng)流場、結(jié)構(gòu)特性和風(fēng)-橋相互作用的不確定性,對大跨度拱橋渦振進(jìn)行可靠性評估。

渦振是自激、限幅振動,但會造成結(jié)構(gòu)疲勞,嚴(yán)重影響行車的舒適性,當(dāng)阻尼更低時,橋梁渦振危害更大;然而現(xiàn)有橋梁渦振損傷研究忽略了阻尼、風(fēng)速等不確定性因素。本文采用概率方法評估阻尼、風(fēng)速等不確定性因素的影響,以秋浦河大橋為例,使用數(shù)值模擬對π型梁懸索橋的渦振進(jìn)行研究,考慮阻尼的隨機性,擬合出不同阻尼與最大渦振振幅的函數(shù)曲線,用對數(shù)正態(tài)分布來擬合橋梁特定月份的概率密度,最后使用設(shè)計點法+重要抽樣蒙特卡洛(Monte Carlo,MC)法[8]、直接抽樣MC法、Latin超立方抽樣MC法進(jìn)行行車舒適度的概率性評估。

1 渦振模擬原理

將截面渦激振動系統(tǒng)簡化為豎向、扭轉(zhuǎn)的二維振動彈簧-質(zhì)量-阻尼系統(tǒng),如圖1所示,基于此簡化模型進(jìn)行結(jié)構(gòu)渦激振動模擬計算。首先固定模型進(jìn)行計算,使模型截面的數(shù)值模擬進(jìn)行充分繞流,然后突然釋放模型,使其在氣動力作用下運動,流體和固體進(jìn)一步相互作用,進(jìn)而產(chǎn)生流固耦合。圖1中:cv、ct分別為豎彎、扭轉(zhuǎn)阻尼;kv、kt分別為豎彎、扭轉(zhuǎn)剛度。

圖1 π型梁橋截面二維振動簡化模型

本文采用交錯方法[9]進(jìn)行流固耦合,此方法在傳統(tǒng)方法基礎(chǔ)上添加流場力的預(yù)測和修正2個步驟,免去迭代,且靈活性高,可以組合各種流體和固體求解器,本文使用數(shù)值模擬平臺“小牛CFD”[10],采用有限元法求解結(jié)構(gòu),有限體積法求解流體。在1個時間步里,該方法包括以下4個計算步驟。

(1)

其中,Δtn、Δtn-1為時間步。

(2)

其中,β為模型參數(shù)。當(dāng)結(jié)構(gòu)與流體密度比很大時,β取值為0.5;當(dāng)結(jié)構(gòu)與流體密度比較小時,β取值范圍為0<β≤0.5。

數(shù)值模擬選用剪切應(yīng)力輸運(shear-stress transport,SST)k-ω(k為湍流動能,ω為特征耗散速度)湍流模型,采用HOPE(High-order Potential Energy)動網(wǎng)格技術(shù)[11],對流固耦合進(jìn)行求解,識別結(jié)構(gòu)渦振響應(yīng)。計算區(qū)域外邊界為矩形,流場參數(shù)設(shè)置如下:流場左側(cè)為速度入口條件,指定速度大小;流場右側(cè)為出口邊界條件,不指定速度(模擬自然流出的邊界);上、下兩側(cè)為對稱邊界條件,其法向速度值等于0,不指定切向速度(模擬自然流出的邊界)。網(wǎng)格劃分方法:主梁截面周邊及尾流附近采用密網(wǎng)格,計算域周邊采用離散網(wǎng)格,中間區(qū)域進(jìn)行合理過渡。網(wǎng)格采用三角形非結(jié)構(gòu)化網(wǎng)格,滿足近壁面網(wǎng)格尺寸Y+<7的要求。

2 π型梁渦激振動模擬與分析

2.1 有限元建模及結(jié)構(gòu)動力特性分析

在ANSYS平臺建立秋浦河大橋三維有限元模型,主梁、橋塔采用兩節(jié)點空間梁單元模擬,主纜、吊桿采用兩節(jié)點三維桿單元模擬,橋面鋪裝、欄桿等荷載采用Mass21單元模擬,主梁橫截面示意如圖2所示(單位為m)。

圖2 秋浦河大橋π型主梁橫截面示意

橋梁結(jié)構(gòu)自振特性分析是進(jìn)行結(jié)構(gòu)動力分析的基礎(chǔ),是檢驗結(jié)構(gòu)模型正確與否的關(guān)鍵環(huán)節(jié)。本文采用子空間迭代法,獲得該橋前50階模態(tài)頻率和振型,受篇幅限制,僅給出前5階模態(tài),見表1所列。

表1 秋浦河大橋前5階模態(tài)

2.2 求解參數(shù)設(shè)置

進(jìn)行網(wǎng)格的無關(guān)性研究,最小網(wǎng)格尺寸為9.833×10-4m,對應(yīng)網(wǎng)格數(shù)量為88 782,并以此為基數(shù),算出各網(wǎng)格數(shù)量與基數(shù)的比值K,計算0°風(fēng)攻角下截面的三分力系數(shù)CD、CL、CM,結(jié)果如圖3所示。由圖3可知,在K>0.4時,三分力系數(shù)變化幅度都在5%以內(nèi),在保證計算精度和效率的前提下,確定K=1.0最為合理。

圖3 三分力系數(shù)隨網(wǎng)格變化曲線

設(shè)定計算風(fēng)速區(qū)間為0.25～3.50 m/s,對應(yīng)該橋?qū)嶋H風(fēng)速為2～28 m/s。在計算風(fēng)速為0.25、0.50、1.00、1.50、1.75、2.00、2.25、2.50、3.00、3.50 m/s下,進(jìn)行豎向渦激振動模擬計算,計算參數(shù)見表2所列。

表2 秋浦河大橋渦激振動模擬計算參數(shù)

2.00 m/s計算風(fēng)速下,Re=6.66×104,網(wǎng)格形式及Y+分布如圖4所示;三角形網(wǎng)格數(shù)為88 782,節(jié)點數(shù)為44 974,采用4層多重網(wǎng)格,最小網(wǎng)格尺寸為9.833×10-4m,邊界條件如圖5所示。

圖4 流場網(wǎng)格和Y+的分布

圖5 計算區(qū)域邊界條件

本文π型主梁高跨比為10,開口率為84.7%,有較低的豎向渦振起振風(fēng)速和較長的豎向渦振鎖定區(qū)間,豎向渦振性能較差[1];文獻(xiàn)[3]研究發(fā)現(xiàn),高跨比8.442、開口率85.1%的π型梁(與本文主梁截面相似)在試驗過程中存在明顯的豎向渦激共振現(xiàn)象。因此,本文主要針對π型主梁豎向渦激振動機理進(jìn)行探究。

2.3 結(jié)果分析

不同計算風(fēng)速v下,無量綱豎向振幅ymax/D(ymax為豎向位移最大值)變化曲線如圖6所示。圖6中，f為渦脫頻率。為了數(shù)據(jù)的可視化,對各個計算風(fēng)速的豎向位移曲線進(jìn)行傅里葉變換,變換后的幅頻曲線如圖7所示(部分風(fēng)速幅頻曲線的振幅過小,振幅(A)需擴大一定倍數(shù))。

圖6 渦振豎向位移隨風(fēng)速變化曲線

圖7 不同計算風(fēng)速下幅頻曲線

由圖6、圖7可知:當(dāng)計算風(fēng)速低于1.00 m/s時,結(jié)構(gòu)位移很小;計算風(fēng)速在1.00～2.50 m/s區(qū)間時,隨著計算風(fēng)速增大,結(jié)構(gòu)位移先快速增大然后減小;在計算風(fēng)速等于1.75 m/s時,豎向位移達(dá)到最大值,為1.16 mm;計算風(fēng)速超過2.50 m/s時,結(jié)構(gòu)位移很小。此現(xiàn)象說明,本文模擬結(jié)果中含有1個渦激鎖定區(qū)間,鎖定計算風(fēng)速在1.00～2.50 m/s。

旋渦脫落頻率整體上隨著風(fēng)速增加而增加,在風(fēng)速鎖定區(qū)間內(nèi),旋渦脫落頻率基本保持在2.201 Hz附近,說明此時結(jié)構(gòu)自振頻率捕獲了旋渦脫落頻率。

0.50 m/s的計算風(fēng)速下,存在渦旋脫落頻率2.202 Hz,但其豎向位移相比于同頻率計算風(fēng)速下的豎向位移小得多。為了解釋上述現(xiàn)象,引入動力放大系數(shù)μ=ymax/yst(yst為靜荷載下的位移),采用結(jié)構(gòu)動力學(xué)里的共振放大系數(shù)進(jìn)行計算,表達(dá)式為:

(3)

其中:fn為結(jié)構(gòu)固有頻率;ξ為阻尼比。

在0.50、1.75 m/s計算風(fēng)速下,對升力進(jìn)行傅里葉變換后的功率譜如圖8所示。由圖8可知,達(dá)到豎彎振動頻率時,0.50 m/s時的功率遠(yuǎn)小于1.75 m/s時的功率。0.50 m/s計算風(fēng)速下的豎向位移較低,推測主要原因為計算風(fēng)速低,輸入的激勵較小,導(dǎo)致主梁本身的靜位移小。

圖8 2種計算風(fēng)速下的功率譜曲線

根據(jù)文獻(xiàn)[1],π型主梁高跨比10、開口率84.7%橋梁起振的折算風(fēng)速(v/(fD))應(yīng)為10,對應(yīng)渦振區(qū)間長度應(yīng)為14～18;本文折算風(fēng)速為9.415,對應(yīng)渦振區(qū)間長度為14.123,與文獻(xiàn)[1]十分接近,驗證了本文數(shù)值模擬的可靠性。

3 橋梁渦振易損性曲線和概率評價

3.1 橋梁運營期阻尼與渦振的關(guān)系

為探究橋梁阻尼比隨運營時間的變化規(guī)律,本文采用汕頭海灣大橋、錢塘江三橋、黃石長江大橋3座橋梁的多年監(jiān)測數(shù)據(jù)[12],繪制3座橋梁阻尼比ξ隨時間變化曲線,如圖9所示。

圖9 3座橋梁阻尼比隨時間衰減過程曲線

阻尼的存在讓原先的物理振動系統(tǒng)隨時間增加逐漸趨于靜止,使系統(tǒng)能量逐漸消散,文獻(xiàn)[13-14]研究發(fā)現(xiàn),不同加載、不同測點等因素影響下,阻尼變化較大。從圖9可以看出,隨時間增長,ξ呈先上升后下降、最后趨于穩(wěn)定的趨勢。阻尼比的改變會對橋梁渦振響應(yīng)產(chǎn)生一定影響,因此本文對不同阻尼比下結(jié)構(gòu)的渦激振動響應(yīng)進(jìn)行研究。

3.2 ymax隨ξv變化規(guī)律

改變結(jié)構(gòu)豎彎阻尼比ξv,其余參數(shù)不變,計算風(fēng)速區(qū)間取0.50～3.00 m/s,對應(yīng)秋浦河大橋?qū)嶋H風(fēng)速u為4～24 m/s,計算步長為0.50 m/s，將計算結(jié)果換算為該橋?qū)嶋H風(fēng)速-振幅,如圖10所示。

圖10 秋浦河大橋?qū)嶋H豎向渦振振幅變化曲線

ymax隨ξv變化計算參數(shù)見表3所列。

表3 ymax隨ξv變化計算參數(shù)

從圖10可以看出，在4種ξv下,π型截面結(jié)構(gòu)均發(fā)生渦激振動現(xiàn)象,實際發(fā)生渦激振動風(fēng)速區(qū)間為8～20 m/s,在u=14 m/s時,結(jié)構(gòu)渦激振動振幅最大。根據(jù)文獻(xiàn)[15],計算主梁豎向基頻對應(yīng)的渦激振動容許振幅為109.215 mm。由圖10、表3可知:當(dāng)ξv=0.004時,結(jié)構(gòu)最大位移超過規(guī)范允許值,渦激振動振幅過大;當(dāng)ξv為規(guī)范建議值0.010時,結(jié)構(gòu)最大位移為55.7 mm;豎向渦振幅值隨阻尼比增大而減小。

為探究豎彎阻尼比ξv與豎向渦振最大幅值ymax的關(guān)系,使用MATLAB中BP神經(jīng)網(wǎng)絡(luò)和擬合函數(shù)命令,其中BP神經(jīng)網(wǎng)絡(luò)的完整訓(xùn)練需要先傳入數(shù)據(jù)給input,通過隱含層的映射在輸出層得到BP神經(jīng)網(wǎng)絡(luò)仿真值,將訓(xùn)練得到的仿真值與目標(biāo)值進(jìn)行對比和誤差分析,修改權(quán)值和閾值,使仿真值不斷逼近目標(biāo)值,對最終結(jié)果進(jìn)行擬合,得到表達(dá)式為:

ymax=aeb ξv+ced ξv

(4)

其中:a=575.30;b=-374.80;c=48.16;d=-15.49。

對(4)式兩邊取對數(shù),使用最小二乘法再次擬合,得到ξv和ymax的關(guān)系式為:

lnymax=a2lnξv+b2

(5)

其中,a2、b2為待定系數(shù)。

根據(jù)(5)式對表3中的lnymax和lnξv進(jìn)行最小二乘法擬合,得到擬合關(guān)系式為:

lnymax=-1.167lnξv-1.343

(6)

由(6)式可知,lnymax和lnξv呈線性關(guān)系,利用(6)式可以推算不同阻尼比下豎向渦激振動的最大振幅,lnymax與lnξv呈負(fù)相關(guān),增大ξv可以減小ymax,可得ξv=5.67×10-3時,實際結(jié)構(gòu)位移達(dá)到規(guī)范允許幅值[15]。為驗證該結(jié)論的準(zhǔn)確性,在v=1.75 m/s、ξv=5.67×10-3下進(jìn)行計算流體動力學(xué)(computational fluid dynamics,CFD)模擬,結(jié)果如圖11所示。由圖11可知,模型振幅為1.78 mm,對應(yīng)實際豎向振幅為106.8 mm,與允許振幅109.215 mm十分接近。

圖11 CFD模擬位移時程曲線

3.3 ξv的概率密度分布

橋梁的阻尼在1年中受氣候等影響，出現(xiàn)很大程度變化。根據(jù)文獻(xiàn)[12]中西堠門大橋連續(xù)6 a的監(jiān)測數(shù)據(jù),得到該橋不同月下的豎彎阻尼比ξv變化情況，如圖12所示。

由圖12可知,在每年的8—12月,ξv變化較小,是具有代表性、偏保守的較低阻尼比。為研究橋梁阻尼的演變規(guī)律,取ξv分布為對數(shù)正態(tài)分布;對每年8—12月的lnξv進(jìn)行擬合,得到其概率分布的Q-Q圖,如圖13所示;ξv概率密度函數(shù)分布如圖14所示。

圖12 西堠門大橋ξv隨時間演變情況

圖13 西堠門大橋ln ξv的Q-Q圖

圖14 西堠門大橋ξv概率密度函數(shù)分布

從圖13可以看出,lnξv的期望正態(tài)值總是圍繞其參照線呈均勻分布,基本符合一條直線,由此可知，西堠門大橋每年8—12月的ξv總是服從對數(shù)正態(tài)分布。本文π型梁橋與西堠門大橋在地理上位于相近緯度,氣候變化較為類似,因此可將西堠門大橋阻尼比方差應(yīng)用于本文π型梁橋。

3.4 概率損傷評價

結(jié)構(gòu)的損傷概率定義為:

(7)

其中:f(ξv)為影響結(jié)構(gòu)性能的隨機變量ξv的概率密度函數(shù);G(ξv)為結(jié)構(gòu)功能函數(shù)。極限狀態(tài)面G(ξv)=0將ξv的設(shè)計空間分為失效域G(ξv)≤0和安全域G(ξv)>0。

根據(jù)(5)式、(7)式,定義功能函數(shù)G(ξv)=z0-(a2lnξv+b2),對于pf較小的點,通過采用設(shè)計點法+重要抽樣MC法來估計結(jié)構(gòu)的損傷概率,其余點使用直接抽樣MC法估計結(jié)構(gòu)的損傷概率;其中,通過設(shè)計點法可找出結(jié)構(gòu)最可能失效的點,通過重要抽樣MC法確定選用正態(tài)分布概率密度函數(shù)抽樣。

(8)

表4 Sperling穩(wěn)定性評價標(biāo)準(zhǔn)

采用重要抽樣MC法模擬1×106次,可得pf=0.05的誤差有99%的概率小于0.364%,而使用直接抽樣MC法,需要模擬9.5×106次。設(shè)N1為直接抽樣MC模擬次數(shù)，N2重要抽樣MC模擬次數(shù)，求解最優(yōu)化問題min(N1-N2),有

(9)

其中:β=Φ-1(1-pf);Vp=ε/μα/2,ε為pf的誤差,μα/2為標(biāo)準(zhǔn)正態(tài)分布的上α/2分位點,1-α為小于ε的置信度。

當(dāng)α=1%,ε=0.364%時,對pf>0.500 5的點,重要抽樣MC法模擬計算量大幅增加,用直接抽樣MC法效果更好。

本文π型梁橋最大振幅下ξv的損傷概率如圖15所示。

從圖15可以看出：隨著ξv不斷減小,結(jié)構(gòu)的損傷概率不斷增加;當(dāng)ξv=7.78×10-3,橋梁渦振有5%的概率超過W=2.00下的允許振幅,結(jié)構(gòu)具有95%的可靠度;當(dāng)ξv=8.00×10-3時,不同W取值下的損傷概率均小于5%,滿足結(jié)構(gòu)渦振抗風(fēng)要求。

圖15 本文π型梁橋最大振幅下ξv的損傷概率變化曲線

以上損傷概率不考慮當(dāng)?shù)貥蛭伙L(fēng)速,可視為渦振最大幅值的損傷概率;下面結(jié)合當(dāng)?shù)貥蛭伙L(fēng)速發(fā)生概率,分析渦振區(qū)間的損傷概率。

匯總1979—2020年安徽望江地區(qū)逐日最大風(fēng)速,得到對應(yīng)擬合的非參數(shù)化概率密度,如圖16所示。根據(jù)功能函數(shù)14±6[exp(a2lnξv+b2)-exp(z0)]/exp(a2lnξv+b2),使用Latin超立方抽樣MC法,計算在不同W、ξv下幅值超過評價標(biāo)準(zhǔn)的風(fēng)速區(qū)間,當(dāng)W=2.00時,ymax隨ξv、u的變化規(guī)律如圖17所示。

圖16 望江地區(qū)逐日最大風(fēng)速概率密度

由圖17可知,超過評價標(biāo)準(zhǔn)的鎖定風(fēng)速區(qū)間隨ξv增大而減小。

圖17 ymax隨ξv、u變化規(guī)律

結(jié)合擬合風(fēng)速的概率密度,把風(fēng)速區(qū)間代入累積概率密度中,得到逐年ξv的損傷概率,即考慮渦振區(qū)間的損傷概率pa,如圖18所示。

圖18 渦振區(qū)間下ξv的損傷概率變化曲線

基于表4中行車舒適度的評價標(biāo)準(zhǔn),本文將結(jié)構(gòu)損傷程度分為2個等級:等級Ⅰ,W=2.00對應(yīng)損傷概率為5%的ξv;等級Ⅱ,W=2.50對應(yīng)損傷概率為5%的ξv。由圖15、圖18可知,與只考慮渦振最大幅值損傷概率pf相比，渦振區(qū)間下對應(yīng)損傷概率pa均有不同程度下降,結(jié)果更加保守。

綜上所述,ξv在每年8—12月服從對數(shù)正態(tài)分布,逐年有整體下降趨勢;對于本文π型梁橋,建議每年8—12月測定ξv,當(dāng)ξv降到等級Ⅰ時,加強監(jiān)測頻率,當(dāng)ξv進(jìn)一步降到等級Ⅱ時,加強養(yǎng)護(hù)頻率,添加人工阻尼。

4 結(jié) 論

本文基于大變形動網(wǎng)格和有限體積法實現(xiàn)π型截面的流固耦合模擬,可以準(zhǔn)確識別渦激共振的發(fā)生區(qū)間和共振頻率;在此基礎(chǔ)上,研究π型梁渦激共振下的損傷評價。主要結(jié)論如下:

(1) 通過π型主梁截面的渦激振動數(shù)值模擬可知,渦脫頻率整體上隨著風(fēng)速增加而增加,渦脫頻率在2.201 Hz被結(jié)構(gòu)俘獲;豎向渦激振幅在鎖定區(qū)間內(nèi)先增大后減小。

(2) π型梁懸索橋的豎向渦激共振現(xiàn)象在相當(dāng)廣泛的阻尼比范圍內(nèi)均存在,隨著阻尼比減小,結(jié)構(gòu)豎向渦激振動振幅最大值逐漸增大,但豎向渦振鎖定風(fēng)速區(qū)間不變;增大阻尼比可以有效減小豎向渦激振幅。

(3) 受氣候等影響,阻尼在1年中有較大變化,每年8—12月豎彎阻尼比服從對數(shù)正態(tài)分布;基于行車舒適性要求和阻尼比變化,使用設(shè)計點法+重要抽樣MC法、直接抽樣MC法,成功擬合了最大幅值的損傷概率曲線;依據(jù)當(dāng)?shù)貥蛭伙L(fēng)速,使用Latin超立方抽樣MC法,擬合渦振區(qū)間的損傷概率曲線,最后提出基于行車舒適度的渦振2級損傷評價方法。