基于貝葉斯估計(jì)的振動基頻提取及預(yù)測算法*

王代君 陳 星 龔惠云 熊 瑋

(1.江西省交通運(yùn)輸科學(xué)研究院有限公司 南昌 330000； 2.江西省交通監(jiān)控指揮中心 南昌 330000)

拉索受力狀態(tài)是斜拉橋結(jié)構(gòu)健康評估的重要依據(jù)[1-2]。當(dāng)前使用較多的索力測量方法有油壓表讀數(shù)法、壓力傳感器法、磁通量法、振動頻率法等[3]，其中磁通量法和振動頻率法均可適用于動態(tài)在線監(jiān)測。振動頻率法無論是在理論體系還是傳感器開發(fā)上都已經(jīng)比較成熟，成為當(dāng)前應(yīng)用最為廣泛的索力測量技術(shù)[4]。然而，振動頻率法要依賴于振動的基頻，而不同的拉索其工況各不相同，信號的質(zhì)量也參差不齊，因此基頻的準(zhǔn)確提取成為實(shí)際應(yīng)用過程中的主要難點(diǎn)。以往需要在后期對振動模態(tài)進(jìn)行人工分析才能準(zhǔn)確地得到振動基頻，關(guān)于基頻在線提取的算法研究并不多見。

貝葉斯方法在工程應(yīng)用中有很強(qiáng)的實(shí)用性，近年來，國內(nèi)外許多學(xué)者將其應(yīng)用于橋梁結(jié)構(gòu)分析并取得了良好的應(yīng)用效果。溫雨桐[5]、董俊等[6]、秦超等[7]、Xie Y.N.等[8]、樓燦洪等[9]將貝葉斯方法應(yīng)用在橋梁模態(tài)參數(shù)識別和可靠性分析等領(lǐng)域，李源等[10]研究了基于貝葉斯更新的斜拉橋施工過程索力預(yù)測，該方法針對的是施工過程，未涉及運(yùn)營期的長期監(jiān)測。本文的創(chuàng)新點(diǎn)在于，在斜拉橋運(yùn)營期長期監(jiān)測的應(yīng)用場景下，利用貝葉斯方法預(yù)測未來一段時期的索力變化趨勢，以期為斜拉橋的養(yǎng)護(hù)和應(yīng)急管理提供參考。

1 振動頻率法基本原理

振動頻率法是基于弦振動理論提出的，利用索力與自振頻率之間的關(guān)系來推導(dǎo)索力。根據(jù)弦振動理論的觀點(diǎn)，斜拉索的自由振動方程可表示為

(1)

式中：x為沿索長方向的坐標(biāo)；t為時間坐標(biāo);v(x,t)為t時刻索上各點(diǎn)的豎向位移；EI為索的抗彎剛度；T為索力；m=W/g為索的線密度，W為單位索長的重力，g為重力加速度。

如果斜拉索兩端采用鉸支，可解得索力為

(2)

式中：fn為索的第n階自振頻率；l為索長。

在不考慮抗彎剛度的情況下，索力解可進(jìn)一步簡化為

(3)

式中：K=4Wl2/g為比例系數(shù)?？梢?，在拉索的材料和生產(chǎn)工藝一致的條件下，K與l索長成正比。對于拉索的基頻(1階自振頻率)，有

T=K(f1)2

(4)

不難發(fā)現(xiàn)，求解索力的關(guān)鍵是要先求出基頻f1。聯(lián)合式(3)和式(4)，有

(5)

(6)

根據(jù)式(5)和式(6)可以發(fā)現(xiàn)索振動的2個重要特性。

特性一。在同一根索的同一次索力測量中，拉索振動功率譜上的各個峰值是等間距的，且間距的寬度值與基頻相等。

特性二。在同一根索的同一次索力測量中，拉索振動功率譜上的所有高階頻率都是基頻的整數(shù)倍，且第n階自振頻率等于基頻的n倍。

2 振動基頻提取

索力計(jì)算的關(guān)鍵在于基頻的提取，這也是當(dāng)前工程領(lǐng)域的難點(diǎn)之一。在實(shí)際應(yīng)用中，由于受到噪聲的影響，頻譜曲線中通常存在大量的雜亂頻率成份，索的自振頻率甚至?xí)辉肼曀蜎]。

根據(jù)前述索振動的兩重要特性均可實(shí)現(xiàn)振動基頻的間接提取。對于特性一，只要將拉索振動功率譜上的所有高階頻率提取出來，然后分析各頻率的差值即可。對于特性二，在提取出高階頻率之后，通過分析某一階頻率屬于第幾階來達(dá)到基頻識別的目的。但某階頻率所屬的階數(shù)往往是未知的，因此從實(shí)際工程應(yīng)用的角度來看，第一個特性的可操作性更強(qiáng)。

基頻提取算法過程示意圖見圖1。

圖1 基頻提取算法過程示意圖

首先獲取1組振動原始時域信號，見圖1a)，然后對時域信號進(jìn)行快速傅里葉變換，得到振動功率譜見圖1b)?？梢悦黠@看到，功率譜中的低階自振頻率和高階自振頻率均被噪聲所淹沒，只有中頻段幾組較強(qiáng)的自振頻率仍比較明顯。將這幾組較強(qiáng)的自振頻率提取出來，見圖1c)。將后一組自振頻率與前一組進(jìn)行作差，可知其頻率差是相等的，均為0.937 5 Hz。

以上基頻提取算法主要利用了拉索的第一個振動特性，本文將其稱為“差頻法”。差頻法要求傳感器有較好的激勵效果，且原始信號具有較高的信噪比，否則算法的性能會明顯下降，仍然需要人工介入才能求解出第n階自振頻率或相鄰兩階的頻率差。例如，在以下情況下將無法采用差頻法求解索力。

1) 如果只有1個自振頻率被成功提取出來，則“作差”就無從談起，直接導(dǎo)致基頻提取失敗。而這種情況在工程中是很常見的。

2) 當(dāng)2階自振頻率之間存在較大的干擾頻率時，該干擾頻率會被誤判為某一階自振頻率，使得頻差不再是一個固定的常數(shù)，致基頻求解失敗。

為解決這一缺陷，本文采用樸素貝葉斯的方法對算法進(jìn)行優(yōu)化。

3 貝葉斯優(yōu)化

3.1 貝葉斯定理

貝葉斯定理可表示為

(7)

式中：π為參數(shù)的概率分布；π(θ)為先驗(yàn)概率；π(θ|x)為后驗(yàn)概率；f(x|θ)為觀測樣本分布，稱為似然函數(shù)(Likelihood)，對先驗(yàn)概率起到修正作用；積分項(xiàng)稱為標(biāo)準(zhǔn)化常量，積分區(qū)間Θ為參數(shù)θ的值域。

3.2 參數(shù)模型分析

為了應(yīng)用貝葉斯定理對基頻求解算法進(jìn)行優(yōu)化，首先需要明確基頻的先驗(yàn)概率分布。在拉索的全生命周期中，其基頻是一個連續(xù)變化的物理量，而在實(shí)際測量中，只能測得一個離散的基頻序列，為了還原其概率分布，對已測得的歷史基頻值進(jìn)行直方圖統(tǒng)計(jì)，發(fā)現(xiàn)其滿足高斯分布，實(shí)測基頻值主要在其理論值附近波動，離理論值越遠(yuǎn)，其概率就越小。因而本文采用正態(tài)分布N(μ,σ2)對其進(jìn)行擬合，從而得到連續(xù)的先驗(yàn)概率密度曲線。

根據(jù)貝葉斯定理可知，單純依靠先驗(yàn)概率對事件進(jìn)行預(yù)測是不準(zhǔn)確的，還需要通過似然函數(shù)對先驗(yàn)概率進(jìn)行修正。本文采用最近1 d的基頻值作為觀測樣本求解似然函數(shù)。實(shí)驗(yàn)表明，基頻值在短周期內(nèi)仍然服從高斯分布，因而也可以采用正態(tài)分布模型N(μ′,σ′2)求解似然函數(shù)。

為了提高算法的魯棒性，在采用差頻法的同時，還運(yùn)用了“主頻法”求解基頻。所謂“主頻法”就是利用拉索振動的第二個特性，找出拉索振動功率譜上能量最大的高階自振頻率(即主振頻率fmax)，然后再除以其階數(shù)n，即可求得基頻。利用該方法，即使只提取到1個自振頻率，也可以較可靠地求出基頻，彌補(bǔ)了差頻法的不足。

3.3 貝葉斯決策

在正態(tài)分布模型中，離中心點(diǎn)越近的樣本出現(xiàn)的概率越高，如果差頻法和主頻法計(jì)算的結(jié)果不一致，取離中心點(diǎn)最近的結(jié)果作為基頻測量值即可。根據(jù)正態(tài)分布的性質(zhì)，下一次樣本出現(xiàn)在期望附近的概率最高。顯然，要使式(7)取得最大值，只要使f(x|θ)π(θ)取得最大值極可。即下一次的基頻預(yù)測值為

(8)

3.4 模型的進(jìn)化

在計(jì)算出當(dāng)前基頻值后，將其真實(shí)值作為新的“歷史”數(shù)據(jù)對先驗(yàn)概率分布模型和似然函數(shù)進(jìn)行更新，重新確定模型參數(shù)，通過動態(tài)的學(xué)習(xí)使模型始終保持最優(yōu)參數(shù)，不斷提高基頻提取及預(yù)測的準(zhǔn)確性。

基于貝葉斯估計(jì)的算法優(yōu)化流程圖見圖2。

圖2 基于貝葉斯的算法優(yōu)化流程圖

4 工程應(yīng)用案例

4.1 工程概況

藥都大橋位于江西省宜春市樟樹市，全橋長2 529 m，寬29.2 m，其中主橋長(168 m+400 m+168 m)736 m。藥都大橋?yàn)殡p塔雙索面鋼混疊合梁斜拉橋，主塔為寶瓶形，塔高124.52 m，全橋拉索共84對168根空間扇形雙索面斜拉索，采用扭絞型平行鋼束拉索，主跨400 m，為贛江主跨最大橋梁。為滿足橋梁長期自動監(jiān)測的需要，藥都大橋配備了橋梁結(jié)構(gòu)健康監(jiān)測系統(tǒng)。系統(tǒng)采用單向加速度傳感器結(jié)合頻譜法對拉索索力情況進(jìn)行監(jiān)測。傳感器安裝于橋梁各拉索位置3.5 m處，全橋索力測點(diǎn)共56個。傳感器采集的原始物理量為電壓，采集頻率為25 Hz，經(jīng)數(shù)據(jù)處理后計(jì)算出加速度和索力值。

4.2 數(shù)據(jù)分析

1) 數(shù)據(jù)建模與迭代。以主跨上游15號索CA02-15為例，對該索2020年1-2月份的振動基頻測量值進(jìn)行直方圖統(tǒng)計(jì)并進(jìn)行高斯函數(shù)建模，得到基頻的先驗(yàn)概率分布為N(0.950 25,0.002 8282)，見圖3a)。對2020年2月23-29日的振動基頻測量值進(jìn)行直方圖統(tǒng)計(jì)并進(jìn)行高斯函數(shù)建模，得到基頻的似然函數(shù)模型為N(0.951 92,0.002 2942)，見圖3b)。

在初始模型的基礎(chǔ)上，使用新求解出的基頻值對模型進(jìn)行更新，如圖3a)所示的虛線為經(jīng)過243次迭代后的概率密度曲線。在模型迭代的過程中，參數(shù)也會不斷自我調(diào)整，見圖4。

2) 性能分析。根據(jù)模型對基頻進(jìn)行提取并預(yù)測，所得結(jié)果見圖5a)，預(yù)測誤差見圖5b)。不難看出，預(yù)測曲線與測量曲線的變化趨勢是一致的，且總體預(yù)測誤差在±1.5%以內(nèi)。

圖3 主跨上游15號索數(shù)據(jù)建模

圖4 先驗(yàn)概率模型參數(shù)迭代

圖5 預(yù)測值與測量值對比及預(yù)測誤差

為了評價(jià)模型的誤差情況，可以采用采用均方根誤差RMSE、平均絕對值誤差MAE、平均絕對百分比誤差MAPE等指標(biāo)進(jìn)行描述，其中RMSE既可表征絕對誤差，又可以表征誤差的波動性，因此本文采用RMSE對模型的誤差和穩(wěn)定性進(jìn)行分析。RMSE定義如式(9)。

(9)

將誤差曲線分為[1,80]、[81,160]、[161,243]3個區(qū)間，分段計(jì)算誤差曲線的RMSE，結(jié)果分別為0.371 3、0.323 3、0.179 9，說明模型經(jīng)過訓(xùn)練后穩(wěn)定性不斷提高，驗(yàn)證了模型的有效性。

4.3 全橋基頻提取

全橋共56個索力測點(diǎn)，基頻提取結(jié)果見圖6a)。根據(jù)式(4)計(jì)算全橋索力并與成橋索力進(jìn)行對比，得到索力偏差見圖6b)。

根據(jù)圖6a)顯示的全橋基頻分布可以發(fā)現(xiàn)，越靠近橋塔的拉索，其基頻越大，反之亦然。這是因?yàn)榭拷鼧蛩睦鏖L度較短，遠(yuǎn)離橋塔的拉索長度更長，根據(jù)索長與基頻成反比的特性可知，計(jì)算的基頻分布與工程實(shí)際完全相符。從圖6b)顯示的結(jié)果不難看出，測量值與成橋索力對比的偏差均在±10%以內(nèi)，根據(jù)T/CECS 529-2018《大跨度橋梁結(jié)構(gòu)健康監(jiān)測系統(tǒng)預(yù)警閾值標(biāo)準(zhǔn)》，全橋索力均處于正常狀態(tài)。

圖6 全橋基頻及索力偏差

5 結(jié)語

本文對基于貝葉斯估計(jì)的振動基頻提取及預(yù)測算法進(jìn)行了研究，結(jié)論如下。

1) 基于弦振動理論推導(dǎo)了索力與振動基頻及各級自振頻率之間的關(guān)系，證明在同一根索的同一次索力測量中，振動功率譜上的各階自振頻率是等間距的，且間距的寬度值與基頻相等。所有高階頻率都是基頻的整數(shù)倍，且第n階自振頻率等于基頻的n倍。

2) 分別利用拉索的2個重要振動特性設(shè)計(jì)了“差頻法”和“主頻法”，采用2種方法聯(lián)合求解基頻，彌補(bǔ)了差頻法的不足，這是本文基頻提取算法最大的優(yōu)勢。

3) 采用高斯函數(shù)對先驗(yàn)概率分布和似然函數(shù)進(jìn)行建模，根據(jù)貝葉斯定理對基頻進(jìn)行預(yù)測，充分利用歷史數(shù)據(jù)蘊(yùn)含的有效信息提升預(yù)測模型的精度，使總體預(yù)測誤差控制在±1.5%以內(nèi)。經(jīng)過迭代后的模型參數(shù)能自適應(yīng)調(diào)整，并逐漸趨于穩(wěn)定，驗(yàn)證了模型的有效性。這也是本文的主要創(chuàng)新點(diǎn)。

4) 根據(jù)基頻計(jì)算的全橋索力分布與工程實(shí)際相符，與成橋索力對比的總體偏差在±10%以內(nèi)，全橋索力均處于正常狀態(tài)。

5) 貝葉斯方法是一種基于概率論的預(yù)測方法，需要有足夠多的先驗(yàn)信息才能達(dá)到較高的準(zhǔn)確性，因此本文提供的方法對于長期索力預(yù)測具有一定的優(yōu)勢。