王 超， 毛 羚

(1.湖北工業(yè)大學(xué)土木建筑與環(huán)境學(xué)院 武漢， 430034) (2. 武漢理工大學(xué)土木工程與建筑學(xué)院 武漢，430070)

引 言

對于土木工程結(jié)構(gòu)損傷診斷、健康監(jiān)測與安全評估來說，模態(tài)參數(shù)識別是一個非常重要的研究課題。最近幾十年，各國學(xué)者提出了許多模態(tài)參數(shù)識別技術(shù)，如頻域分解法、時域分解法、隨機(jī)子空間識別、基于響應(yīng)傳遞比的參數(shù)識別法[1]和基于貝葉斯理論的模態(tài)參數(shù)識別法[2]等。這些參數(shù)識別方法都是針對時不變系統(tǒng)。然而，實(shí)際土木工程結(jié)構(gòu)在運(yùn)營過程中由于環(huán)境或者荷載變化可能導(dǎo)致結(jié)構(gòu)發(fā)生時變，如不同時刻構(gòu)件不均勻溫度變換會引起結(jié)構(gòu)模態(tài)頻率隨時間改變；橋上車輛荷載的變化引起結(jié)構(gòu)拉索索力的變化，導(dǎo)致拉索剛度時變；車輛過橋時由于車輛位置不斷改變引起車橋系統(tǒng)質(zhì)量時變。因此許多學(xué)者提出了時變參數(shù)識別方法進(jìn)行系統(tǒng)識別，如狀態(tài)空間小波識別法[3]、S變換法[4]、HHT(Hilbert-Huang transform, 簡稱HHT)法[5]、自適應(yīng)稀疏時頻分析法[6]、同步擠壓小波變換法[7]和解析模態(tài)分解法(analytical modal decomposition簡稱AMD)[8]等。

變分模態(tài)分解(variational mode decomposition,簡稱VMD)是最近提出的一種新的信號自適應(yīng)分解方法[9],與HHT類似，它可以將非平穩(wěn)的時變信號自適應(yīng)的分解成多個具有有限帶寬的模態(tài)分量信號之和，但是它沒有HHT分解時遇到的模態(tài)混疊問題。由VMD分解得到模態(tài)分量需再采取合適方法提取信號的瞬時頻率。HHT是采用Hilbert變換的方法進(jìn)行分析得到Hilbert譜，但是它會產(chǎn)生負(fù)頻率和端點(diǎn)效應(yīng)問題。利用小波變換技術(shù)也可以提取信號瞬時頻率，通常采用的復(fù)Morlet小波并不是完全解析的，當(dāng)中心頻率較小時可能產(chǎn)生負(fù)頻率泄露，導(dǎo)致時頻分析結(jié)果產(chǎn)生偏差。

因此，筆者采用VMD和廣義Morse小波(generalized Morse wavelets，簡稱GMW)相結(jié)合的技術(shù)來識別時變結(jié)構(gòu)瞬時頻率。首先,利用信號在小波時頻域具有稀疏性的特性，引入快速閾值迭代技術(shù)對信號進(jìn)行去噪處理;然后,采用VMD將去噪信號分解為多個模態(tài)分量;最后,采用GMW對每個模態(tài)分量進(jìn)行時頻分析提取對應(yīng)的瞬時頻率。通過數(shù)值算例和一個移動小車試驗(yàn)驗(yàn)證了提出的方法。

1 快速閾值迭代降噪

現(xiàn)場實(shí)測信號通常會受到噪聲污染,這里引入快速閾值迭代算法(fast iterative shrinkage-thresholding algorithm, 簡稱FISTA)進(jìn)行降噪。含噪信號通過小波變換可以轉(zhuǎn)化到時間-尺度域，噪聲能量將分布到非常廣泛的時間-尺度空間，而有效信號能量將集中在一個局部的時間-尺度域空間，具有稀疏特性，因此，通過最優(yōu)稀疏性可以實(shí)現(xiàn)將噪聲和信號盡量分離。假設(shè)實(shí)測的含噪信號y表示為

y=x+n=Iw(C)+n

(1)

其中:x為真實(shí)無噪信號;C為x小波變換系數(shù)矩陣;Iw(·)為小波逆變換;n為噪音。

稀疏性對應(yīng)于l0范數(shù)最小化問題，對于非凸復(fù)雜性優(yōu)化問題，可以采用l1范數(shù)來替代求解[10],因此，上述問題可以轉(zhuǎn)換為求解優(yōu)化問題

(2)

其中：s.t.為約束條件;ε為小量。

式(2)有約束優(yōu)化問題可以通過施加拉格朗日乘子正則化方法轉(zhuǎn)變?yōu)闊o約束優(yōu)化問題來求解[11]

(3)

其中:λ為拉格朗日乘子。

這里采用Beck等[12]提出的快速閾值迭代算法(IST)通過下式快速迭代求解式(3)的解

(4)

其中：WT(·)為小波變換;Tλ(·)為軟閾值函數(shù);k為迭代數(shù);t為迭代參數(shù)。

其中：c為小波系數(shù)矩陣C中的小波系數(shù)。

為了提高收斂效率，λ按指數(shù)衰減變化取值

(7)

其中：N為迭代求解次數(shù);λmax=0.99max|C0|;λmin=max|C0|/q，q為系數(shù)C0的數(shù)量。

通過上述算法對響應(yīng)數(shù)據(jù)進(jìn)行分析，可以有效降低噪音干擾。

2 變分模態(tài)分解

VMD是一種新的信號自適應(yīng)分解技術(shù)，將實(shí)信號分解為具有最優(yōu)特定稀疏特性的多個分量信號(模態(tài))之和。這里的稀疏特性采用分解信號的頻譜帶寬來衡量， VMD要求分解的每個模態(tài)分量具有圍繞各自中心頻率最緊的帶寬。

為了估計每個模態(tài)分量的帶寬，VMD利用希爾伯特變換計算各模態(tài)對應(yīng)的解析信號, 得到其單邊頻譜, 然后通過頻率平移將每個模態(tài)頻譜平移到各自估計的中心頻率位置，最后計算頻移后信號梯度范數(shù)的平方來估計各模態(tài)分量的帶寬。這樣最優(yōu)模態(tài)分解轉(zhuǎn)化為如下有約束的變分優(yōu)化問題來實(shí)現(xiàn)[13]

(8)

其中：{si}和{ωi}分別為分解的模態(tài)分量和對應(yīng)的中心頻率;δ(t)為狄拉克函數(shù);*為卷積;f為分解的信號。

通過引入二次罰函數(shù)項(xiàng)和拉格朗日乘子ζ來處理上式的重構(gòu)約束條件，拉格朗日乘子ζ用來施加嚴(yán)格的約束條件，而二次罰函數(shù)項(xiàng)可以提高收斂性能，將式(8)的約束優(yōu)化問題轉(zhuǎn)變?yōu)槭?9)無約束的函數(shù)優(yōu)化問題

(9)

其中：α為帶寬約束參數(shù)。

α增加可能導(dǎo)致分解的模態(tài)分量帶寬減小，中心頻率誤差估計增大；α減小可能導(dǎo)致分解的模態(tài)帶寬增加，包含更多的噪音成分。α建議取值區(qū)間為100～10 000, 文中取2 000。

通過交替方向乘子法(alternating direction method of multipliers,簡稱ADMM)可迭代求取原變分優(yōu)化問題的解。步驟如下：

(10)

其中：F-1(·)表示傅里葉逆變換。

(11)

4) 按式(12)進(jìn)行迭代更新

(12)

循環(huán)迭代步驟2～4直到滿足收斂條件

(13)

其中:ζ為迭代收斂容差，文中取1e-7。

3 廣義Morse小波瞬時頻率識別

對于N自由度的時變結(jié)構(gòu)系統(tǒng)，其振動響應(yīng)z(t)可分解為N階振型模態(tài)響應(yīng)的疊加和

(14)

其中:zm(t)為第m階振型模態(tài)響應(yīng)。

每階模態(tài)響應(yīng)可表示為時變幅值A(chǔ)m(t)和相位φm(t)的乘積

xm(t)=Am(t)eiφm(t)

(15)

通過VMD對結(jié)構(gòu)響應(yīng)信號z(t)進(jìn)行分析，可以將結(jié)構(gòu)不同階振型模態(tài)分離成各模態(tài)分量，然后分別對各模態(tài)分量采用連續(xù)小波變換進(jìn)行時頻分析，提取對應(yīng)的小波脊線，由小波尺度與信號的頻率關(guān)系可以得到結(jié)構(gòu)的瞬時頻率[7]。

連續(xù)小波時頻分析通常采用復(fù)Morlet小波，然而，復(fù)Morlet小波只是當(dāng)中心頻率較大時近似解析，當(dāng)中心頻率較小時可能產(chǎn)生負(fù)頻率泄露，導(dǎo)致分析結(jié)果產(chǎn)生一定偏差。文中引入GMW[13]進(jìn)行時頻分析, 該小波是完全解析的，不存在負(fù)頻率泄露問題。GMW頻域定義為

Ψβ, γ(ω)=H(ω)2(eγ/β)β/γ)ωβe-ω γ

(16)

其中：H(ω)為階躍函數(shù);β和γ為控制小波時域和頻域衰減的參數(shù)。

當(dāng)γ=3時，GMW接近高斯函數(shù)，具有很好的時頻分辨率。文中采用的GMW參數(shù)為γ=3和β=27，對應(yīng)的時域和頻域波形如圖1所示(示意圖，幅值無量綱)。

圖1 GMW時域和頻域波形Fig.1 Wave form of GMW in time domainand frequency domain

4 數(shù)值算例

為驗(yàn)證提出的算法, 首先采用如下具有不同時變頻率成分組合成的振動信號進(jìn)行數(shù)值驗(yàn)證

(17)

信號采樣頻率為100 Hz，并添加10%的高斯白噪聲。這里噪聲強(qiáng)度10%定義為噪聲的2-范數(shù)與信號的2-范數(shù)之比。

首先,利用FISTA對含噪信號進(jìn)行降噪處理, 降噪前后信號的時頻分布圖(又稱小波量圖，尺度坐標(biāo)為指數(shù)坐標(biāo)，無量綱)如圖2所示。

圖2 信號小波量圖Fig.2 Signal wavelet scalogram

圖中顏色亮度表示小波系數(shù)幅值大小，越亮表示越大，對應(yīng)能量就越大?？梢钥闯?，信號成分在小波域分布具有稀疏性，比較集中，能量較大，而噪聲分布較廣，通過FISTA處理后噪聲明顯減少。

隨后,采用VMD對降噪后的信號進(jìn)行分析，分解得到如圖3所示的4個模態(tài)分量，原信號各成分也顯示在圖中用于比較。

圖3 VMD分解的各模態(tài)分量(黑線)和對應(yīng)的原信號xi(藍(lán)線)Fig.3 The decomposed mode components by VMD (dark line) and corresponding signal xi (blue line)

從圖可以看出,VMD將原信號不同成分完好的分解開來，其中第2個成分x2具有突變特性，VMD將其突變前后分解成兩個部分，分解的各部分與原信號對應(yīng)各成分吻合的非常好。

采用GMW對VMD分解的各模態(tài)分量進(jìn)行連續(xù)小波分析提取對應(yīng)的小波脊線，識別其時變頻率。另外，也采用HHT的方法的識別結(jié)果用于對比，原信號的瞬時頻率的理論值用來驗(yàn)證識別結(jié)果的準(zhǔn)確性，分析結(jié)果如圖4所示。

圖4 信號各成分時變頻率識別結(jié)果Fig.4 The identified time-varying frequency of each component

由圖4可知，HHT方法識別的3個信號成分瞬時頻率受噪音干擾波動較大，對于信號突變成分突變時刻存在一定的偏差，而且HHT方法的端點(diǎn)效應(yīng)影響很大。提出的方法3個信號成分瞬時頻率識別結(jié)果與理論值吻合的很好，只在端部由于小波端點(diǎn)效應(yīng)稍有偏離，但是遠(yuǎn)小于HHT的影響。

5 移動荷載試驗(yàn)

5.1 試驗(yàn)?zāi)Ｐ?/h3>

為進(jìn)一步驗(yàn)證提出的方法，設(shè)計了一個移動小車試驗(yàn)，采用馬達(dá)牽引小車通過一個簡支鋼板梁來模擬移動荷載作用下的結(jié)構(gòu)時變特性。板梁尺寸為1 000 mm×150 mm×9 mm，質(zhì)量為10.555 kg，采用圓鋼桿支撐模擬的簡支支撐。為避免小車進(jìn)入和離開主梁時產(chǎn)生大的振動與沖擊，在板梁首尾分別連接了一段過渡的引梁和尾梁，但它們與主梁斷開。小車設(shè)有兩軸，軸距為0.117 m，質(zhì)量為2.308 kg，小車內(nèi)通過砝碼調(diào)整重量, 本試驗(yàn)時小車總重為5.296 kg。主梁上設(shè)置了一根與主梁中線平行的穿過小車的的鋼絲，使小車能沿板梁中線移動。試驗(yàn)?zāi)Ｐ腿鐖D5所示。

圖5 試驗(yàn)?zāi)Ｐ蛨DFig. 5 The experimental mode

在梁底板起點(diǎn)和終點(diǎn)，1/4,1/2,3/4跨處布設(shè)了應(yīng)變片，在節(jié)點(diǎn)3,5,7,9處(見圖6)布設(shè)了加速度計，測點(diǎn)布置如圖6所示。

圖6 測點(diǎn)布置圖Fig. 6 The layout of measure points

5.2 模態(tài)試驗(yàn)

首先進(jìn)行模態(tài)試驗(yàn)。在節(jié)點(diǎn)4處用力錘錘擊主梁，采集結(jié)構(gòu)的加速度響應(yīng)，采樣頻率為1 000 Hz。進(jìn)行模態(tài)分析并對結(jié)構(gòu)的有限元模型進(jìn)行校正。隨后將小車當(dāng)做移動質(zhì)量塊，利用校正后的有限元模型計算小車重心通過不同節(jié)點(diǎn)位置時結(jié)構(gòu)的第1階模態(tài)頻率，結(jié)果如表1所示。

表1 小車重心通過不同節(jié)點(diǎn)時結(jié)構(gòu)第1階頻率(FEA)Tab.1 The first frequency of structure when the vehicle cross different node (FEA)

5.3 移動小車試驗(yàn)

通過馬達(dá)牽引小車勻速通過試驗(yàn)梁，小車移動時導(dǎo)致主梁上質(zhì)量分布發(fā)生改變, 引起結(jié)構(gòu)的時變。采集響應(yīng)的應(yīng)變和加速度響應(yīng)，采樣頻率為200 Hz。節(jié)點(diǎn)1和11處測得的應(yīng)變響應(yīng)如圖7所示。

圖7 節(jié)點(diǎn)1和節(jié)點(diǎn)11處測點(diǎn)應(yīng)變響應(yīng)Fig.7 Response of strain by sensor at node 1 and node11

通過首尾節(jié)點(diǎn)應(yīng)變測點(diǎn)響應(yīng)、主梁長度、小車軸距可以計算小車通過主梁時速度。應(yīng)變峰值對應(yīng)小車前后輪正好經(jīng)過測點(diǎn)位置，由圖7測得前輪通過節(jié)點(diǎn)1的時刻為4.925 s，后輪通過節(jié)點(diǎn)11的時刻為7.895 s。由此估算得到小車速度為0.376 m/s。

采用提出的方法對節(jié)點(diǎn)9處加速度測點(diǎn)測得的響應(yīng)進(jìn)行分析識別結(jié)構(gòu)的第一階瞬時頻率。為進(jìn)行比較，也分別采用HHT和復(fù)Morlet連續(xù)小波變換方法(continuous wavelet transform,簡稱CWT)進(jìn)行時變頻率識別，結(jié)果如圖8所示。

圖8 節(jié)點(diǎn)9處測點(diǎn)識別的瞬時頻率Fig.8 Identified instantaneous frequency by sensor at node 9

從識別的結(jié)果可以看出，采用提出了VMD和GMW方法識別的結(jié)果比較光滑，性能較好，通常的CWT識別結(jié)構(gòu)波動稍大，主要存在不同模態(tài)間干擾和噪音干擾，HHT識別結(jié)果最差，端點(diǎn)效應(yīng)影響很大。

采用GMW方法提取小波脊線時，對識別精度的影響主要是信號小波脊線對應(yīng)的各階成分是否被充分激起，另外信號各階成分頻率間距太小，間距過密可能會造成能量的相互干擾。

隨后對所有測點(diǎn)(節(jié)點(diǎn)3,5,7,9)測得的加速度響應(yīng)進(jìn)行分析，首先用FISTA算法進(jìn)行降噪，然后采用VMD和GMW相結(jié)合的方法識別結(jié)構(gòu)的第1階瞬時頻率。由應(yīng)變估計的小車移動速度和通過主梁的起止時間，可以推導(dǎo)出小車在不同時刻在主梁上的位置。結(jié)合表1中有限元計算的結(jié)果，可以得到小車通過主梁不同時刻的第1階瞬時頻率有限元計算結(jié)果，與提出的方法進(jìn)行對比分析，如圖9所示。

圖9 所有測點(diǎn)識別的結(jié)構(gòu)瞬時頻率Fig.9 Identified instantaneous frequency by all sensors

由圖可見，所有測點(diǎn)數(shù)據(jù)的識別結(jié)果與有限元計算結(jié)果吻合的較好，主要在兩端位置稍有偏差，分析偏差的原因:a.尾梁跟主梁連接位置還是存在一定沖擊振動干擾;b.小波端點(diǎn)效應(yīng)有一定影響。另外，所有測點(diǎn)識別結(jié)果相差很小，表明提出的方法可以有效識別結(jié)構(gòu)時變參數(shù)。

6 結(jié)束語

土木工程結(jié)構(gòu)在運(yùn)營過程中通常表現(xiàn)出時變特性，筆者提出了基于VMD和GMW的方法識別時變結(jié)構(gòu)瞬時頻率。對于實(shí)測信號的噪聲污染問題，提出采用快速閾值迭代算法(FISTA)進(jìn)行降噪處理。數(shù)值模擬結(jié)果表明，提出的FISTA算法通過將信號轉(zhuǎn)換到小波域利用信號稀疏特性可以較好地對信號進(jìn)行降噪，VMD和GMW法可以有效準(zhǔn)確地識別時變結(jié)構(gòu)瞬時頻率。移動小車試驗(yàn)進(jìn)一步驗(yàn)證了方法的有效性和精度，提出的方法比HHT和CWT具有更好的效果。