基于貝葉斯壓縮感知的復(fù)合材料結(jié)構(gòu)沖擊載荷識別研究

嚴 剛, 孫 浩

(1.南京航空航天大學(xué)機械結(jié)構(gòu)力學(xué)及控制國家重點實驗室， 江蘇 南京 210016；2. 麻省理工學(xué)院土木與環(huán)境工程系， 美國 馬薩諸塞州 坎布里奇 MA 02139)

引 言

現(xiàn)代飛行器結(jié)構(gòu)大量使用先進復(fù)合材料，但復(fù)合材料結(jié)構(gòu)在低能量外物的沖擊下，容易造成目不可檢的內(nèi)部損傷，導(dǎo)致結(jié)構(gòu)強度及穩(wěn)定性的急劇下降，嚴重威脅飛行安全[1]。因此利用集成在結(jié)構(gòu)中的先進傳感元件，準確監(jiān)測到外物低速沖擊并重建沖擊載荷時間歷程，對及時做出檢測和維護決策，提高復(fù)合材料結(jié)構(gòu)的安全具有重要意義[2]。

在復(fù)合材料結(jié)構(gòu)沖擊載荷時間歷程識別領(lǐng)域已經(jīng)發(fā)展的技術(shù)和方法中，除少數(shù)采用神經(jīng)網(wǎng)絡(luò)等機器學(xué)習(xí)方法對大量樣本進行訓(xùn)練之外[3]，大部分可歸類為基于模型的方法。如Chang等基于假設(shè)模態(tài)法建立了加筋復(fù)合材料結(jié)構(gòu)在沖擊載荷作用下的動力學(xué)模型，通過一種最優(yōu)平滑濾波算法來重建沖擊載荷時間歷程[4]；嚴剛等基于Chang等的模型，結(jié)合參數(shù)化的沖擊載荷時間歷程，采用微種群遺傳算法優(yōu)化識別了參數(shù)，從而近似重建沖擊載荷時間歷程[5]；Hu等通過有限元方法建立了沖擊載荷與結(jié)構(gòu)響應(yīng)之間的關(guān)系，采用二次規(guī)劃最優(yōu)算法獲得沖擊載荷時間歷程[6]。

在基于模型的沖擊載荷時間歷程識別中，利用脈沖響應(yīng)函數(shù)在時域進行反卷積的方法被認為是一類數(shù)學(xué)上直接的方法。但時域反卷積方法的難點在于其不適定性，通常需要進行正則化來獲得穩(wěn)定有界的解[7]。Tikhonov正則化是當(dāng)前在沖擊載荷識別中普遍使用的一種正則化方法[8-11]。在Tikhonov正則化中，一個困難是正則化因子的選取，其選擇是否適當(dāng)直接影響到正則化解的效果[12]。Jin等人在貝葉斯統(tǒng)計推斷的理論框架下，提出了一種增廣的Tikhnov正則化方法(又稱貝葉斯正則化方法)，能根據(jù)測量數(shù)據(jù)自適應(yīng)地確定最優(yōu)正則化因子[13]。Yan和Sun成功地將之應(yīng)用于復(fù)合材料結(jié)構(gòu)沖擊載荷識別，取得了很好的識別效果[14]。

近年來，壓縮感知和稀疏正則化方法在信號和圖像處理領(lǐng)域取得了很大的成功，因此也被研究者應(yīng)用于結(jié)構(gòu)沖擊載荷識別。如Qiao等采用不同的“字典”來稀疏表示沖擊載荷，通過稀疏正則化方法重構(gòu)這些“字典”的系數(shù)來實現(xiàn)對沖擊載荷的重建，表明了利用稀疏正則化對沖擊載荷識別是可行有效的[15]。但傳統(tǒng)的稀疏正則化同Tikhonov正則化類似，也需要事先確定正則化因子。本文針對這一問題，提出一種基于貝葉斯壓縮感知的新方法對沖擊載荷時間歷程進行識別，并通過對復(fù)合材料夾層結(jié)構(gòu)進行的實驗研究來驗證所提出方法的有效性和可應(yīng)用性。

1 沖擊載荷識別模型

在線彈性范圍內(nèi)，結(jié)構(gòu)上任意點的輸入-輸出關(guān)系都可由第一類積分方程表示。假設(shè)只有一個沖擊載荷f(t)作用于復(fù)合材料結(jié)構(gòu)，在結(jié)構(gòu)表面布置了傳感器來測量對應(yīng)的沖擊響應(yīng)ys(t)。若忽略沖擊載荷作用點處局部材料退化等因素，可以將復(fù)合材料結(jié)構(gòu)視為一時不變系統(tǒng)，沖擊載荷與傳感點沖擊響應(yīng)之間的關(guān)系可表示為如下的卷積

(1)

式中hs(t)為沖擊點與傳感點之間的脈沖響應(yīng)函數(shù)。

對式(1)在時域進行離散化處理，可獲得如下的矩陣形式

Ys=HsF

(2)

其中，

(3)

(4)

(5)

F, Ys,Hs分別為時域離散化后的沖擊載荷向量、沖擊響應(yīng)向量和傳遞函數(shù)矩陣，式中f(n)，ys(n)和hs(n)分別是f(t)，ys(t)和hs(t)在離散時刻tn=nΔt處的值，其中Δt是采樣周期，Nsz是采樣點數(shù)。

式(2)是沖擊載荷時間歷程識別的基本方程，但由于存在模型誤差和測量噪聲，式(2)中的傳遞函數(shù)矩陣Hs通常是不適定的，具有很大的條件數(shù)。因此對式(2)直接求逆來識別沖擊載荷，通常得到的解是不穩(wěn)定的，具有非常大的誤差，需要引入正則化來獲得近似的穩(wěn)定有界解。

2 稀疏正則化系數(shù)重構(gòu)問題

為識別沖擊載荷時間歷程，本文首先提出采用一組共M個高斯核函數(shù)作為基函數(shù)來表示未知沖擊載荷。高斯核函數(shù)的形式為

(6)

F=GQ

(7)

將式(7)代入式(2)，則沖擊響應(yīng)向量與基函數(shù)系數(shù)向量之間的關(guān)系可以表示為

Ys=TQ

(8)

其中

Τ=HsG

(9)

對于沖擊載荷，采用高斯核函數(shù)表示后，基函數(shù)系數(shù)向量具有高度稀疏性，因此沖擊載荷時間歷程識別問題可轉(zhuǎn)換為下式描述的基于L2-L1范數(shù)的基函數(shù)系數(shù)稀疏正則化重構(gòu)問題

(10)

3 貝葉斯壓縮感知

本文采用貝葉斯壓縮感知來對式(10)描述的稀疏正則化問題進行求解，其優(yōu)點是可以避免式(10)中正則化因子的選取[16]。采用貝葉斯建模方法，將系數(shù)向量Q封裝在后驗概率密度函數(shù)p(Q|Ys,α,σ2)中，由貝葉斯定理可得[17]

(11)

式中p(Ys|Q,σ2)是似然函數(shù)，p(Q|α)是未知系數(shù)向量的先驗概率密度函數(shù)，p(Ys|α,σ2)是邊緣似然函數(shù)。

假設(shè)測量獲得的沖擊響應(yīng)與模型輸出之間的誤差滿足均值為零、方差為σ2的高斯分布(σ2未知)，似然函數(shù)可寫為

(12)

同樣定義Q的先驗分布為高斯分布，寫為

(13)

將式(12)和(13)代入式(11)，可以得到

(14)

為確定超參數(shù)向量α和誤差方差σ2，采用相關(guān)向量機方法進行求解。這歸結(jié)為一個第二類型的最大似然問題，其形式為[18]

(15)

其中，Ω=σ2I+TΛ-1TT。式(15)的詳細算法可參閱文獻[18]。

將獲得的超參數(shù)向量α和誤差方差σ2代入解向量Q的后驗概率密度函數(shù)p(Q|Ys,α,σ2)，并以p(Q|Ys,α,σ2)最大處對應(yīng)的值作為未知系數(shù)向量Q的最大后驗概率解，代入式(7)就可實現(xiàn)沖擊載荷時間歷程的識別。

4 實驗研究

4.1 實驗裝置與試驗件

為驗證本文所提出方法的有效性和可應(yīng)用性，建立了如圖1所示的實驗系統(tǒng)，在一復(fù)合材料夾層結(jié)構(gòu)上進行實驗研究。該夾層結(jié)構(gòu)呈六邊形，主要尺寸如圖2所示，上下表面面板采用T300/3234碳纖維預(yù)浸料層合而成，上面板厚度為5 mm，下面板厚度3 mm，夾芯層采用PMI泡沫，厚度為17 mm。在上面板表面粘貼了4個直徑為10 mm，厚度為1 mm的PZT壓電片，分別命名為S1, S2, S3, S4，組成傳感網(wǎng)絡(luò)，用于監(jiān)測沖擊載荷作用下結(jié)構(gòu)的動態(tài)應(yīng)變。建立了如圖2所示的直角坐標(biāo)系，各個PZT壓電片的坐標(biāo)在表1中列出。實驗中考慮了6個隨機選定的沖擊點，分別命名為IP1, IP2, IP3,IP4, IP5, IP6，其坐標(biāo)在表2中列出。實驗時將復(fù)合材料夾層結(jié)構(gòu)固定在支架上，采用PCB 086C03力錘在結(jié)構(gòu)上施加沖擊載荷，沖擊載荷信號通過B&W 61001功率放大器供電放大后被NI PXIe-5105數(shù)字化儀采集，結(jié)構(gòu)在沖擊載荷作用下的響應(yīng)由PZT壓電片感應(yīng)后同樣由數(shù)字化儀采集，采樣頻率設(shè)置為50 kHz。

圖1 實驗裝置與試驗件Fig.1 Experimental setup and test structure

圖2 復(fù)合材料夾層結(jié)構(gòu)示意圖Fig.2 Illustration of sandwich composite structure

S1S2S3S4坐標(biāo)/mm(50, 350)(350, 300)(300, 0)(0, 100)

表2 沖擊點位置信息

實驗分為兩個階段，第一階段采用模態(tài)實驗方法獲得各沖擊點與各傳感器之間的脈沖響應(yīng)函數(shù)。首先應(yīng)用力錘在各沖擊點連續(xù)沖擊3次，同時記錄PZT壓電片測量到的動態(tài)響應(yīng)信號，經(jīng)傅里葉變換，獲得沖擊點與各傳感器之間的平均頻響函數(shù)，然后再由逆傅里葉變換獲得傳遞函數(shù)矩陣。第二階段，再次對各沖擊點進行沖擊，記錄PZT壓電片獲得的沖擊響應(yīng)，用于對沖擊載荷時間歷程進行識別。圖3所示為力錘作用在IP1時的沖擊載荷時間歷程及對應(yīng)的沖擊響應(yīng)信號，其余位置的沖擊載荷時間歷程與沖擊響應(yīng)信號與此類似。

圖3 力錘作用在沖擊點IP1時信號Fig.3 Measured signals when impact hammer acts on impact point IP1

4.2 實驗結(jié)果

為了對沖擊載荷時間歷程進行識別，需要確定作為基函數(shù)的高斯核函數(shù)。本文中高斯核函數(shù)的個數(shù)設(shè)為M= 501，每個核函數(shù)寬度設(shè)置為均勻相等，取hk=T/(M-1)，其中T為沖擊響應(yīng)信號的時間長度，本文中取T= 20 ms。每個核函數(shù)的中心位置均勻地分布在整個時間域上，即tk=(k-1)hk(k= 1,2,…,M)。圖4所示為第101個高斯核函數(shù)的示意圖。

圖5 傳感器S1響應(yīng)識別沖擊點IP1處沖擊載荷結(jié)果Fig.5 Identification results for impact point IP1 by S1 response

圖4 基函數(shù)中第101個高斯核函數(shù)示意圖Fig.4 101st Gaussian kernal function

接著對所采集到的各組沖擊響應(yīng)信號進行分析來識別對應(yīng)的沖擊載荷時間歷程。為了比較識別效果，分別采用直接求逆法、Tikhonov正則化和貝葉斯正則化3種方法進行對比。圖5(a)所示為沖擊載荷作用在沖擊點IP1時，采用傳感器S1的響應(yīng)信號直接求逆識別的結(jié)果。從圖中可以看出，由于傳遞函數(shù)矩陣的不適定性，直接求逆的解是不穩(wěn)定的，超出實際沖擊載荷幾個數(shù)量級。圖5(b)所示為采用Tikhonov正則化和貝葉斯正則化對同一響應(yīng)信號進行識別的結(jié)果，其中Tikhonov正則化因子由L曲線法確定[12]。由圖中可以看出，L曲線選擇的正則化參數(shù)雖然消除了解的不穩(wěn)定性，但識別的沖擊載荷低于實際的沖擊載荷，表現(xiàn)為明顯的過正則化；而貝葉斯正則化從測量數(shù)據(jù)中自適應(yīng)地確定了正則化因子，在一定程度上識別出了沖擊載荷的主要部分，但其識別結(jié)果尾部具有很強的振蕩性，且和Tikhonov正則化識別結(jié)果一樣存在負分量，仍需進一步改進。

圖6 本文方法識別IP1處沖擊載荷結(jié)果(采用的響應(yīng)來自S1, S2, S3, S4)Fig.6 Identification results for impact point IP1 with the proposed method(by responses from S1, S2, S3, S4)

圖7 采用不同傳感器響應(yīng)識別沖擊點IP1處沖擊載荷誤差Fig.7 Error of the identification results for impact point IP1 by different sensor responses

圖6(a)所示為沖擊載荷作用在沖擊點IP1時，采用本文所提出方法對傳感器S1的響應(yīng)信號進行識別的結(jié)果。從圖上可以看出，識別的沖擊載荷時間歷程比較平滑，沒有明顯的振蕩，與實際的沖擊載荷吻合得非常好，識別的沖擊載荷的形狀、時間長度以及峰值基本與實際的一致。圖6(b)～(d)所示分別為沖擊載荷作用在沖擊點IP1時，本文所提出方法對其他3個傳感器信號進行識別的結(jié)果，從圖中可以看出，都取得了很好的識別效果。

同樣，對其他沖擊點的沖擊載荷進行識別時，本文所提出方法也都能獲得類似的識別結(jié)果，與實際沖擊載荷很好地吻合。為了進一步比較識別結(jié)果的準確性，定義相對誤差為

(16)

以及峰值誤差為

(17)

式中Ftrue和Freconstructed分別表示實際和識別的沖擊載荷時間歷程。圖7(a)和(b)所示是不同方法分別利用4個傳感器響應(yīng)信號識別沖擊點IP1處沖擊載荷的相對誤差和峰值誤差比較。

圖8(a)和(b)所示是不同識別方法對不同沖擊點處沖擊載荷識別的平均相對誤差和平均峰值誤差比較(4個傳感器響應(yīng)識別誤差的平均值)。從圖8中可以看出，與Tikhonov正則化以及貝葉斯正則化方法相比，本文所提出方法顯著降低了識別結(jié)果的誤差，在本文實驗研究的6個沖擊點沖擊載荷識別結(jié)果中，最大相對誤差僅為12.25%，最大峰值誤差僅為7.23%，表明本文所提出方法識別的沖擊載荷與實際沖擊載荷基本一致，驗證了其有效性和可應(yīng)用型。

圖8 不同沖擊點處識別沖擊載荷的平均誤差Fig.8 Average error of the identification results for different impact points

5 結(jié) 論

本文提出了一種基于貝葉斯壓縮感知的新方法，對作用在復(fù)合材料結(jié)構(gòu)上的沖擊載荷時間歷程進行識別。在采用一組高斯核函數(shù)近似表示未知沖擊載荷時間歷程后，將沖擊載荷重建問題轉(zhuǎn)換為基函數(shù)系數(shù)的稀疏正則化重構(gòu)問題，通過貝葉斯壓縮感知方法來求解該稀疏正則化重構(gòu)問題，獲得基函數(shù)的最優(yōu)稀疏系數(shù)，繼而重建沖擊載荷時間歷程。

對復(fù)合材料夾層結(jié)構(gòu)的實驗研究結(jié)果表明：本文提出的方法能準確地識別沖擊載荷時間歷程，識別的沖擊載荷與實際的沖擊載荷吻合非常好，其形狀、時間長度及峰值與實際沖擊載荷基本一致，定義的相對誤差小于15%，峰值誤差小于10%，均顯著低于Tikhonov正則化和貝葉斯正則化的識別誤差，表明了本文所提出方法的有效性和可應(yīng)用性。

致謝 感謝南京航空航天大學(xué)航空宇航學(xué)院王新峰副教授提供的復(fù)合材料夾層結(jié)構(gòu)用于本文的實驗研究。

[1] 沈 真, 楊勝春, 陳普會.復(fù)合材料抗沖擊性能和結(jié)構(gòu)壓縮設(shè)計許用值 [J]. 航空學(xué)報, 2007, 28(3): 561-566.

Shen Zhen, Yang Shengchun, Chen Puhui. Behaviors of composite materials to withstand impact and structural compressive design allowableness [J].ActaAeronautica Et Astronautica Sinica, 2007, 28(3): 561—566.

[2] Staszewski W J, Mahzan S, Traynor R. Health monitoring of aerospace composite structures - active and passive approach [J]. Composites Science and Technology, 2009, 69(11-12): 1678—1685.

[3] Ghajari M, Sharif-Khodaei Z, Aliabadi M H, et al. Identification of impact force for smart composite stiffened panels [J]. Smart Materials and Structures, 2013, 22(8): 085014.

[4] Seydel R, Chang F K. Impact identification of stiffened composite panels: I. system development [J]. Smart Materials and Structure, 2001,10(1): 354—369.

[5] 嚴 剛, 周 麗.加筋復(fù)合材料結(jié)構(gòu)的沖擊載荷識別 [J]. 航空學(xué)報, 2008, 29(5): 1150—1156.

Yan Gang, Zhou Li. Impact load identification for stiffened composite structure [J]. Acta Aeronautica Et Astronautica Sinica, 2008, 29(5): 1150—1156.

[6] Hu N, Fukunaga H, Matsumoto S, et al. An efficient approach for identifying impact force using embedded piezoelectric sensors [J]. International Journal of Impact Engineering, 2007, 34(7): 1258—1271.

[7] Inoue H, Harrigan J J, Reid S R. Review of inverse analysis for indirect measurement of impact force [J]. Applied Mechanics Review, 2001, 54(6): 503—524.

[8] Jacquelin E, Bennani A, Hamelin P. Force reconstruction: analysis and regularization of a deconvolution problem [J]. Journal of Sound and Vibration, 2003, 265(1): 81—107.

[9] 王 鋒, 武 龍, 張小慶, 等.動態(tài)載荷辨識問題的正則化求解新方法 [J].振動工程學(xué)報, 2016, 29(1): 31—37.

Wang Feng, Wu Long, Zhang Xiaoqing, et al. A technique for solving dynamical force identification problems by Tikhonov regularization method [J]. Journal of Vibration Engineering, 2016, 29(1): 31—37.

[10] 馬 超, 華宏星.一種基于新的正則化技術(shù)的沖擊載荷識別法 [J].振動與沖擊, 2015, 34(12): 164—168.

Ma Chao, Hua Hongxing. Impact force identification based on improved regularization technique [J]. Journal of Vibration and Shock, 34(12): 164—168.

[11] 常曉通, 閆云聚, 劉 鎏, 等.基于Green函數(shù)和正則化的動態(tài)載荷識別方法 [J].振動、測試與診斷, 2014, 34(1): 124—129.

Chang Xiaotong, Yan Yunju, Liu Liu, et al. Applications of regularization method in dynamic load identification [J]. Journal of Vibration, Measurement & Diagnosis, 2014, 34(1): 124—129.

[12] Hansen P C. Analysis of discrete ill-posed problems by means of the L-curve [J]. SIAM review, 1992, 34(4): 561—580.

[13] Jin B, Zou J. A Bayesian inference approach to the ill-posed Cauchy problem of steady-state heat conduction [J]. International Journal of Numerical Methods in Engineering, 2008, 76(4): 521—544.

[14] Yan G, Sun H, Büyük?ztürk O.Impact load identification for composite structures using Bayesian regularization and unscented Kalman filter [J]. Structural Control and Health Monitoring, 2017, 24(5): e1910.

[15] Qiao B J, Zhang X W, Wang C X, et al. Sparse regularization for force identification using dictionaries [J]. Journal of Sound and Vibration, 2016, 368: 71—86.

[16] Ji S, Xue Y, Carin L. Bayesian compressive sensing [J]. IEEE Transactions on Signal Processing, 2008, 56(6): 2346—2356.

[17] 茆詩松, 湯銀才.貝葉斯統(tǒng)計 [M].北京: 中國統(tǒng)計出版社, 2012.

Mao Shisong, Tang Yincai. Bayesian Statistics [M]. Beijing: China Statistics Press, 2012.

[18] Tipping M E. Sparse Bayesian learning and the relevance vector machine [J]. Journal of Machine Learning Research, 2001, 1: 211—244.