基于CEEMDAN-SOBI對橋梁監測撓度的分離研究*

譚冬梅，陳方望，周 強，吳 浩

(1.武漢理工大學 道路橋梁與結構工程湖北省重點實驗室，湖北 武漢 430070；2.華中師范大學 城市與環境科學學院，湖北 武漢 430079)

0 引言

橋梁結構在運營過程中承受多種荷載作用，因此對結構進行撓度監測得到的實時撓度是多項因素共同作用的結果，據已有研究表明，影響因素主要包括車載、溫度作用及混凝土的收縮徐變等[1]。實時監測撓度反映是各項因素共同作用的結果，因此有必要對監測數據進行處理得到各項因素單獨作用下橋梁結構的響應。

針對橋梁結構監測撓度的分離問題，目前已有部分學者對監測得到的響應時程數據，提出不同的方法：針對活載效應的分離，梁宗保等[2]研究小波多尺度分解在活載分離的應用。針對溫度效應的分離，陳夏春等[3]提出應用多元線性回歸模型提取橋梁結構中的溫度效應；劉綱等[4-6]分別提出基于粒子群優化算法的自適應濾波方法、最小二乘支持向量機以及多最小二乘支持向量機來獲取監測數據的溫度效應；孫雅瓊等[7-8]通過現場試驗確定橋梁結構動應變與溫度的關系，提出采用時變多元線性擬合的方法得到溫度效應時程曲線；Tang等[9-11]通過傳統的經驗模態分解(EMD)對信號進行升維，然后進行獨立分量分析；譚冬梅等[12-13]分別采用改進PCA算法結合集成經驗模態分解(EEMD)和MEEMD對實測撓度的日溫差效應進行分離。針對各撓度成分的分離，劉夏平等[14]研究了基于奇異值分解的盲信號分離方法；陳國良等[15]結合時間小波分解和時間序列分析中的中心移動平均法分解得出各成分撓度。

本文針對各撓度成分的分離，提出1種基于CEEMDAN-SOBI的撓度成分分離算法。首先利用CEEMDAN將單通道多源的撓度信號分解為一系列IMF分量，減輕了傳統EMD和EEMD在分解時的模態混疊問題；然后利用排列熵算法對分解得到的各分量進行排列熵值計算，將排列熵值接近的分量重組得到新的分量，采用基于K-L散度的虛假分量判別法識別并剔除重組分量中虛假重組分量，避免直接基于K-L散度剔除虛假模態容易產生的過度剔除。

1 基本理論

1.1 CEEMDAN算法

CEEMDAN算法由TORRES[16]等提出，相比EMD和EEMD算法[17]，該算法更好抑制了模態混疊的問題，可消除IMF分量中的噪聲，自適應獲取集成平均次數和噪聲幅值，提高計算速度[18]。CEEMDAN算法步驟如下：

1)對于第1階IMF分量，其求取方法與EEMD分解方法一致，對原始信號x(t)添加不同的噪聲I次，利用EMD分解I次并求取平均值得到IMF1(t)：

(1)

2)計算第1個殘余量信號R1(t)：

R1(t)=x(t)-IMF1(t)

(2)

3)對信號R1(t)+ε1Ei(ni(t))進行N次EMD分解，得到第2個IMF分量：

(3)

4)對于k∈[2,K]，依次計算第k個殘余信號Rk(t)：

Rk(t)=Rk-1(t)-IMFk(t)

(4)

5)對于每一個Rk(t)，重復步驟3)的過程，得到第k+1個IMF分量為：

(5)

6)重復上述步驟，當殘差信號不能再被繼續分解時，得到最終殘差R(t)：

(6)

由此可得各階IMF分量IMFk(t)和R(t)。

1.2 排列熵

排列熵是用來判斷時間序列復雜程度的1種方式，相比其他常用方法，該方法計算效率高、對噪聲的魯棒性更好[19-20]。序列的排列熵計算步驟如下：

1)已知原始序列z(n)(n=1,2,…,N)，設定維數為m，時間延時為τ，嵌入維度m的取值一般為3～7，嵌入時延τ的取值一般為1。將z(n)構造成如下式矩陣：

(7)

式中：M=N-(m-1)τ，M是矩陣的行數。

2)將矩陣中每1行視為1個向量zj={z(j),z(j+τ),…,z(j+(m-1)τ)},其中j=1,2,…,M，將向量中元素值的大小按照升序排列如下：

z(j+(j1-1)τ)≤…≤z(j+(jm-1)τ)

(8)

式中：j1,j2,…,jm為向量zj中每個元素在排序前的索引號，對于向量中元素值大小相等的情況，即若jp

3)對于每一個向量zj都對應1個符號向量sj={j1,j2,…,jm}，設1個m維重構相空間對應的符號序列的概率密度分別為P1,P2,…,Pk，則對于一維時間序列z(n)的M個重構向量對應的符號序列，排列熵可表示為：

(9)

序列的排列熵值表征了序列的復雜程度，排列熵值越大，則序列越復雜。

1.3 K-L散度

在概率論和信息論中K-L散度[21]也稱之為相對熵，其是對2個概率分布P與Q之間關系的描述，定義為：

(10)

(11)

(12)

式中：P(y)為密度函數的核密度估計；K(·)為核函數；h為給定的正數，通常稱為窗寬或平滑參數。在實際情況下，信號的分布是單峰和對稱的，由上式可得到對稱形式下P和Q之間的K-L散度的定義：

D(P,Q)=DKL(P||Q)+DKL(Q||P)

(13)

在虛假IMF分量識別過程中，K-L散度表征了單通道混合信號與各分量之間關系的緊密程度，K-L散度越小，則分量與原混合信號的關系越緊密，該分量為信號的真實成分，反之關系越疏遠，該分量為需剔除的虛假分量。

1.4 SOBI算法原理

1)計算白化矩陣W，白化矩陣W可根據觀測信號的自相關矩陣的特征值分解得到：

(14)

式中：λ1,λ2,…,λn為原始數據零均值化后的自相關矩陣的n個最大特征值；h1,h2,…,hn分別為特征值對應的特征向量；σ2是噪聲的方差。

2)白化后的信號為：

Z(t)=WX(t)

(15)

3)計算觀測信號白化后的采樣協方差矩陣R(τ)：

R(τ)=E[Z(t+τ)ZT(t)]=ARZ(τ)AT

(16)

τ∈{τj|j=1,2,…,k}

(17)

4)計算正交矩陣U，對于所有的R(τj)，采用聯合近似對角化算法，得出正交矩陣U滿足：

UTR(τj)U=Dj

(18)

式中：Dj是1組對角矩陣。

5)計算分離矩陣和分離信號，混合矩陣A=W-1U，所以分離矩陣B=A-1=UTW，分離信號：

Y(t)=UTWX(t)

(19)

1.5 撓度信號分離流程

撓度分離的步驟：

2)計算每個分量的排列熵，將排列熵值相近的分量相加，組成新的一系列分量；

2 橋梁模擬撓度信號分離

2.1 撓度信號模擬

橋梁在運營過程中的溫度效應的撓度信號主要由日溫差及年溫差組成，橋梁運營過程中還會由于預應力損失及混凝土收縮徐變而引起長期撓度。為更好模擬橋梁撓度信號，在Midas Civil軟件上建立武漢某斜拉橋有限元模型并對結構仿真分析，模型如圖1所示，模擬橋梁撓度信號可表示為：

圖1 武漢某斜拉橋模型Fig.1 Model of a cable-stayed bridge in Wuhan

f(t)=f1(t)+f2(t)+f3(t)

(20)

式中：f(t)為總體撓度信號,mm；f1(t)為日溫差效應,mm，f1(t)=f11(t)+f12(t)，f11(t)整體溫差效應,mm，f12(t)為截面日溫差效應,mm；f2(t)為年溫差效應,mm；f3(t)為長期撓度,mm。

對于溫度作用下撓度信號的模擬，橋梁模型整體升溫1 ℃，跨中向下產生1.86 mm偏移；整體降溫1 ℃，跨中向上產生1.86 mm偏移；主梁截面線性增加溫差1 ℃，跨中向下產生0.54 mm偏移。假定日溫差和年溫差均為正弦變化的周期性函數，其中整體日溫差和截面日溫差共同組成日溫差，取每日的整體日溫差10 ℃，截面日溫差5 ℃以及年溫差30 ℃。由此可得，整體日溫差f11(t)=9.3sin(πt/24)，截面日溫差f12(t)=1.35sin(πt/24)，年溫差f2(t)=27.9sin(πt/8 760)。 對于長期撓度的模擬，參照《公路鋼筋混凝土及預應力混凝土橋涵設計規范》(JTG 3362—2018)，用指數型函數擬合得到長期撓度[13]。各成分撓度信號及總體橋梁撓度模擬信號如圖2所示，為顯示方便，日溫差選取前5 000 h內數據繪圖。

圖2 各成分模擬信號及總撓度模擬信號Fig.2 Simulation signals of each component and total deflection

2.2 撓度信號分離

依據擬定分解流程，對原始信號通過CEEMDAN分解得到13個不同尺度特征的模態分量IMF和1個殘余項R，CEEMDAN算法在一定程度上減輕了EEMD的模態混疊問題，但仍存在一定的模態混疊問題，故需進行進一步處理[15-16]。首先計算各分量的排列熵，將排列熵值相近的分量相加得到重組分量，再結合K-L散度剔除虛假的重組分量。各分量排列熵值如圖3所示，由圖3可知，各分量排列熵值隨著分量頻率的減小而逐漸降低，即頻率低的分量復雜度越低，以各分量排列熵值大小為標準對各分量進行重組。圖3中，IMF1和IMF2的排列熵值明顯大于其他所有分量，且熵值相近，故將兩者相加作為1個新的分量CF1；IMF3，IMF4，IMF5的PE值相差為0.005左右，復雜度相似，故可將三者相加得到新的分量CF2；同理將IMF6，IMF7，IMF8合并為CF3；將IMF9，IMF10，IMF11，IMF12，IMF13合并為CF4；殘余項R的PE值最小，作為重組分量CF5。

圖3 各分量的排列熵值Fig.3 Permutation entropy of each component

對于分解重組得到的分量中必然存在部分虛假分量，本文基于K-L散度理論剔除虛假分量。首先計算各分量與原始信號的K-L散度值，結果如圖4所示，由圖4可知，重組分量CF2，CF4，CF5的K-L散度值最小，故選擇CF2，CF4，CF5作為原始信號的主要成分，其他的予以剔除；最后基于SOBI原理，對提取出來的CF2，CF4，CF5組合為高維矩陣分析，估計原始信號，得到分離結果如圖5所示，為顯示方便，日溫差分離結果選取前5 000 h內數據繪圖。

圖4 基于K-L散度的虛假模態識別Fig.4 False modal recognition based on K-L divergence

圖5 基于SOBI原理的信號分離結果Fig.5 Signal separation results based on SOBI principle

對于分離結果的效果評價，采用分離結果與原始信號間的相關系數及平均絕對誤差來評價，兩者分別如下定義：

相關系數：

(21)

平均絕對誤差：

(22)

表1 各成分分離結果的評價指標Table 1 Evaluation indexes for separation results of each component

3 實測橋梁撓度信號成分分離

3.1 工程背景

武漢某斜拉橋在運營中結合北斗定位導航技術，對其主體結構的三維變形進行實時監測。北斗監測系統較傳統監測方式，具有精度較高、不受氣候條件影響且實時高效率等特點，橋面測點布置如圖6所示。

圖6 武漢某斜拉橋主橋監測測點布置Fig.6 Layout of monitoring points on main bridge of a cable-stayed bridge in Wuhan

考慮相鄰測點布置處的橋面受到溫度作用的影響基本相近，兩相鄰測點記錄的撓度數據中溫度效應及長期撓度的變化趨勢應具有較強的相關性。故對于實測數據的分離效果評價可通過對稱測點分離結果間的相關性來評價。

3.2 實測撓度信號的各成分分離

為驗證本文方法在實際應用中的有效性，采集背景橋梁實時監測撓度1年內的數據，利用CEEMDAN-SOBI對橋梁的各撓度成分進行分離。采集測點編號為BD12和BD35，采集時間段為2017年9月1日00：00：00至2018年8月1日00：00：00，采樣頻率為1 h/次，采集得到的撓度監測數據如圖7所示。為避免部分隨機脈沖噪聲對后續分解的影響，先對原始信號進行高斯平滑，去除異常點，再對信號進行分離。日溫差效應分離結果如圖8所示；年溫差效應分離結果如圖9所示；長期撓度分離結果如圖10所示。

圖7 動撓度實測數據Fig.7 Measured data of dynamic deflection

圖8 日溫差效應分離結果Fig.8 Separation results of daily temperature difference effect

圖9 年溫差效應分離結果Fig.9 Separation results of annual temperature difference effect

圖10 長期撓度分離結果Fig.10 Separation results of long-term deflection

從圖8～10可知，所選跨中對稱測點BD12和BD35在1年內日溫差效應、年溫差效應及長期撓度變化趨勢具有很明顯的相似性，分別計算基于CEEMDAN-SOBI與改進PCA和EEMD[12]方法得到的2測點各撓度成分的分離結果之間的相關系數，如表2所示。從表2可知，各成分間相關系數均較高，證明本文所提方法能較好分離得到實際工程中監測撓度的各成分，且計算精度要優于改進PCA和EEMD方法。

表2 測點BD12和BD35的實測撓度各成分分離結果間的相關系數Table 2 Correlation coefficient between separation results of each component for measured deflection at measuring points BD12 and BD35

4 結論

1)CEEMDAN減輕了傳統EMD和EEMD算法分解信號產生的模態混疊問題，采用排列熵結合K-L散度剔除虛假分量減輕了直接基于K-L散度剔除虛假分量容易產生的過度剔除問題。

2)將CEEMDAN-SOBI算法應用到模擬信號和實測信號的分離，模擬信號分離結果顯示分離值與實際值的相關性及平均絕對誤差趨于理想值；主跨跨中對稱測點實測信號各成分分離結果具有較高的相關性，均在0.94以上，驗證了CEEMDAN-SOBI算法在分離實際橋梁監測撓度各成分的有效性；對比分析了本文計算方法和已有方法對撓度信號的分離效果，結果表明本文方法具有更高精度。

3)本文所提方法是以各效應值時間尺度上的差異為基本條件的，橋梁結構實際運營中，車輛堵車條件下緩慢通行以及積雪等形成的具有與溫度變化尺度相近的荷載會對溫度撓度的準確分離產生一定的干擾，有必要對之監控并剔除。本文提出的方法還可以應用于橋梁長期監測系統中索力和應變的各效應值的分離。