周 玙， 劉 莉， 周思達(dá)， 王大宇， 邵玉佩

(北京理工大學(xué) 宇航學(xué)院 飛行器動力學(xué)與控制教育部重點(diǎn)實(shí)驗(yàn)室，北京 100081)

飛行器結(jié)構(gòu)在飛行或運(yùn)輸過程中經(jīng)常會受到瞬態(tài)載荷的作用。對實(shí)際工作過程中所受載荷進(jìn)行識別，將為飛行器結(jié)構(gòu)的優(yōu)化、故障評估、壽命預(yù)測提供基礎(chǔ)。然而，在實(shí)際工程中，往往載荷無法通過傳感器直接測得，如：測點(diǎn)不可達(dá)、載荷作用點(diǎn)較多，以及涉及流固耦合時(shí)。而一般情況下，結(jié)構(gòu)的動響應(yīng)是較為容易得到的。因此，通過載荷識別方法確定載荷具有較高的實(shí)用價(jià)值。

載荷識別屬于結(jié)構(gòu)動力學(xué)的第二類反問題。在識別瞬態(tài)載荷時(shí)時(shí)域法具有較大優(yōu)勢，國內(nèi)外許多專家學(xué)者對此做過許多研究：周晚林[1]通過有限元逆分析和人工神經(jīng)網(wǎng)路相結(jié)合的綜合方法對智能結(jié)構(gòu)的沖擊載荷進(jìn)行識別；蔡元奇[2]通過假設(shè)載荷形式為二次函數(shù)利用Duhamel積分式推導(dǎo)出了時(shí)域非遞推載荷識別模型；嚴(yán)剛等[3]應(yīng)用遺傳算法對復(fù)合材料沖擊載荷的識別進(jìn)行了研究。劉杰[4]通過求解核函數(shù)矩陣對載荷進(jìn)行識別并研究了其正則化方法；毛伯永等[5]提出了基于瞬態(tài)統(tǒng)計(jì)能量分析理論的沖擊載荷識別方法，識別除了載荷作用位置和輸入能量；張磊等[6]在總體最小二乘法基礎(chǔ) 進(jìn)行Tikhonov正則化后利用共軛梯度法求最優(yōu)解對載荷進(jìn)行識別；張方等[7]利用廣義正交多項(xiàng)式作為形函數(shù)推導(dǎo)了基于時(shí)間有限元的動載荷識別模型，該方法尤其適用于瞬態(tài)載荷的識別。Xu等[8]提出一種根據(jù)不完全測量信息的自適應(yīng)迭代最小二乘法對載荷和參數(shù)進(jìn)行識別。朱斯巖等[9]通過頻域修正時(shí)域載荷識別結(jié)果，并應(yīng)用于實(shí)際工程。

目前，載荷識別主要基于位移響應(yīng)和加速度響應(yīng)。然而，對于運(yùn)行中的飛行器而言，由于傳感器安裝、排線等約束，用于振動測量的加速度計(jì)十分有限，且多普勒激光測振儀(速度)幾乎不可能，因此給其載荷辨識帶來很大困難。

最近，光纖光柵應(yīng)用的興起，提供了輕量、串聯(lián)復(fù)用、抗干擾的結(jié)構(gòu)動應(yīng)變響應(yīng)測試途徑式。然而，應(yīng)變響應(yīng)和結(jié)構(gòu)載荷之間關(guān)系并不直接，需通過應(yīng)變與位移間的關(guān)系間接建立。應(yīng)變到位移是一個空間積分過程，在實(shí)際過程中測點(diǎn)數(shù)量有限，積分結(jié)果誤差很大。目前結(jié)構(gòu)位移模態(tài)和應(yīng)變模態(tài)的識別方法都有一些成熟的方法，許多學(xué)者也推導(dǎo)了應(yīng)變模態(tài)和位移模態(tài)的關(guān)系[10-12]。本文主要圍繞基于應(yīng)變響應(yīng)的沖載荷識別進(jìn)行研究，綜合利用結(jié)構(gòu)的位移模態(tài)和應(yīng)變模態(tài)參數(shù)實(shí)現(xiàn)應(yīng)變到位移的轉(zhuǎn)換避免空間積分，穩(wěn)定、精確地通過應(yīng)變反求出了位移響應(yīng)，進(jìn)而實(shí)現(xiàn)了基于應(yīng)變響應(yīng)的識別。

1 基于位移響應(yīng)的載荷識別

描述多自由度振動系統(tǒng)的方程為

(1)

令Φn×m為系統(tǒng)質(zhì)量歸一化后的振型矩陣，同時(shí)假設(shè)阻尼為比例阻尼。

根據(jù)振型疊加原理

xn×1=Φn×mqm×1

(2)

式中：q為模態(tài)振型坐標(biāo)，代表各振型在振動系統(tǒng)的貢獻(xiàn)度。

將式(2)代入式(1)得

(3)

(4)

得到第j振型坐標(biāo)下的振動方程

(5)

(6)

將式(6)代入式(2)可得在時(shí)刻的動響應(yīng)

(7)

進(jìn)一步可得

(8)

令

式(8)可改寫為

(9)

由于振型坐標(biāo)下的響應(yīng)無法直接獲得，需要將物理坐標(biāo)下的位移響應(yīng)轉(zhuǎn)換到振型坐標(biāo)下，根據(jù)式(2)可得

(10)

(11)

將式(10)代入式(9)可得

(12)

整理得

(13)

(14)

速度信息仍由位移響應(yīng)采用五點(diǎn)差分得到

(15)

2 基于應(yīng)變響應(yīng)的載荷識別

應(yīng)變響應(yīng)相較于位移響應(yīng)對載荷更加敏感，在識別瞬態(tài)載荷時(shí)識別精度更高。同時(shí)應(yīng)變的測量也比位移測量簡便。若利用應(yīng)變響應(yīng)進(jìn)行載荷識別，需建立載荷到應(yīng)變的關(guān)系。一般來說，兩者通過位移相聯(lián)系，應(yīng)變到位移是一個積分過程，實(shí)際工程中測點(diǎn)有限，無法通過應(yīng)變積分得到位移。

本文的基本思路是：通過試驗(yàn)或有限元數(shù)值仿真得到結(jié)構(gòu)的位移模態(tài)和應(yīng)變模態(tài)，找到應(yīng)變和位移的關(guān)系，實(shí)現(xiàn)應(yīng)變到位移的轉(zhuǎn)換。

設(shè)彈性結(jié)構(gòu)的變形位移分別為u，v和w，根據(jù)模態(tài)疊加原理得

(16)

(17)

(18)

式中：qr為時(shí)間和頻率的函數(shù)，它代表的是各階模態(tài)在結(jié)構(gòu)變形中所占的比例，與坐標(biāo)無關(guān)，因此?qr/?x=0。當(dāng)然，根據(jù)上述推導(dǎo)方式，還可以推導(dǎo)出εy，εz以及γxy，γyz，γxz等表達(dá)式。

在三維空間中，設(shè)位移向量為x=[uvw]T，則位移模態(tài)為

(19)

根據(jù)彈性力學(xué)基本原理，應(yīng)變與位移之間的關(guān)系為

(20)

式中：εx，εy和εz為正應(yīng)變，剪應(yīng)變?yōu)棣脁y=αyx+αxy，γyz=αzy+αyz，γzx=αxz+αzx。

根據(jù)模態(tài)疊加原理，得到應(yīng)變張量響應(yīng)表達(dá)式

(21)

考慮到利用應(yīng)變計(jì)或FBG傳感器測量結(jié)構(gòu)應(yīng)變時(shí)，通常只能測量結(jié)構(gòu)正應(yīng)變，因此，只考慮正應(yīng)變時(shí)，則有

(22)

通過比較式(22)和式(19)可以發(fā)現(xiàn)，位移響應(yīng)和應(yīng)變響應(yīng)擁有相同的模態(tài)坐標(biāo)qr，因此，對于同一階模態(tài)對應(yīng)的振型，位移響應(yīng)與應(yīng)變響應(yīng)具有相同的貢獻(xiàn)度。

位移響應(yīng)和應(yīng)變響應(yīng)擁有相同的模態(tài)坐標(biāo)，因此有

(23)

(24)

(25)

式中：ψ和Φ分別為結(jié)構(gòu)的應(yīng)變模態(tài)和位移模態(tài)；Qε和Qx分別為應(yīng)變模態(tài)坐標(biāo)和位移模態(tài)坐標(biāo)，可由最小二乘法解得模態(tài)坐標(biāo)

(26)

由于Qε=Qx，所以得出

(ψTψ)-1ψT·ε=(ΦTΦ)-1ΦT·x

令Sε=(ψTψ)-1ψT，Sx=(ψTΦ)-1ΦT

Sεε=Sxx

(27)

根據(jù)最小二乘法可得應(yīng)變到位移的轉(zhuǎn)換關(guān)系

x=ΦSεε

(28)

為了緩解求逆過程時(shí)的矩陣病態(tài)，對方程進(jìn)行正則化，引入正則化矩陣β，令Sε,β=(ψTψ+β)-1ψT得

x=ΦSε,βε

(29)

當(dāng)以應(yīng)變響應(yīng)作為識別依據(jù)的時(shí)候，根據(jù)位移應(yīng)變關(guān)系得

(30)

值得注意的是在應(yīng)變-位移轉(zhuǎn)換過程中，必須保證應(yīng)變測點(diǎn)nε,d≥m，從而保證式(28)在正定或超定條件下求解。

3 正則化方法及參數(shù)選取

一般在對式(14)和式(28)進(jìn)行求解的時(shí)候采用傳統(tǒng)的最小二乘法就可直接求解。但是當(dāng)系統(tǒng)矩陣條件數(shù)太大時(shí)，系統(tǒng)矩陣病態(tài)，會導(dǎo)致反求結(jié)果失真。此時(shí)需要采用正則化方法來緩解模型的不適定性。

本文參考文獻(xiàn)[13]采用的是擴(kuò)展Tikhonov方法。以一個一般的線性系統(tǒng)為例

AX=Y

(31)

(32)

(33)

正則化參數(shù)λ是平衡正則解的穩(wěn)定性和精度的重要因素，λ過大精度低，過小則不穩(wěn)定。

確定正則化參數(shù)在噪聲水平未知的情況下一般有L曲線法和廣義交叉驗(yàn)證準(zhǔn)則，但是在實(shí)際應(yīng)用時(shí)，L曲線常出現(xiàn)曲線不明顯，參數(shù)不好選擇的情況。所以本文采用廣義交叉驗(yàn)證準(zhǔn)則選取。其形式為

(34)

值得注意的一點(diǎn)就是若其中A是一維向量時(shí)，則采用廣義交叉驗(yàn)證準(zhǔn)則確定的正則化參數(shù)始終為0，即此時(shí)正則化方法無效，此時(shí)需要濾波等手段提高系統(tǒng)的魯棒性。

本文在對式(14)和式(28)進(jìn)行正則化時(shí)均采用上述方法。

4 海洋衛(wèi)星的沖擊載荷識別

為了對所述方法的正確性性和數(shù)值特性進(jìn)行論證，本文根據(jù)文獻(xiàn)[14]建立了海洋一號(HY-1)衛(wèi)星的有限元模型，并嘗試對它進(jìn)行載荷辨識。

模型結(jié)構(gòu)材料都采用碳纖維蜂窩復(fù)合材料(除了對接環(huán))，各板之間采用剛性連接，電子設(shè)備采用非結(jié)構(gòu)質(zhì)量模擬，太陽能帆板采用集中質(zhì)量模擬。

如圖1所示在衛(wèi)星的側(cè)板施加單點(diǎn)瞬態(tài)載荷，載荷形式為

F=2 000sin(π/0.002t),t≤0.02 s

采用有限元軟件Nastran得到衛(wèi)星的位移模態(tài)和應(yīng)變模態(tài)以及動力學(xué)響應(yīng)，測量的位移為各節(jié)點(diǎn)上的位移，測量的應(yīng)變?yōu)楦鲉卧行摹?/p>

圖1 海洋衛(wèi)星有限元模型Fig.1 Oceanic satellite finite element model

4.1 基于位移響應(yīng)的識別

工況描述：模擬衛(wèi)星在軌時(shí)受沖擊情況，衛(wèi)星受到一個垂直于板面的鋸齒波載荷作用(單點(diǎn)激勵)，載荷作用時(shí)間0.04 s，模型不受約束。

單點(diǎn)激勵時(shí)，式(14)中的FI只有一個奇異值，矩陣條件數(shù)為1，矩陣呈良態(tài)，無需正則化，令α=0。

從圖2可以看出，在測量響應(yīng)中無噪聲時(shí)，所建立的載荷識別模型能夠很好地對所作用的沖擊載荷進(jìn)行識別。

從圖3可以看出絕對誤差在開始階段會出現(xiàn)一個極大值，然后迅速穩(wěn)定在某一水平上。開始時(shí)誤差產(chǎn)生極大值的原因是：采用五點(diǎn)差分格式導(dǎo)致速度從第三點(diǎn)開始計(jì)算，而前兩點(diǎn)是通過直接線性插值，速度誤差較大。由于是非遞推關(guān)系，所以此誤差不存在積累，識別結(jié)果的精度在起始幾點(diǎn)后仍得到保證。誤差方向相反大小相等是由于識別模型假設(shè)時(shí)間步長內(nèi)載荷不變，所以一直存在一個與斜率相關(guān)的誤差，載荷斜率方向改變所以誤差方向也隨之改變，用更高階的多項(xiàng)式來近似時(shí)間步長內(nèi)的載荷可以減小該誤差。

圖2 無噪聲識別結(jié)果Fig.2 No noise recognition results

圖3 無噪聲時(shí)識別絕對誤差Fig.3 Absolute error of no noise recognition

4.2 基于應(yīng)變響應(yīng)的識別

整個模型共4 786個節(jié)點(diǎn)，選取150個單元中心的應(yīng)變值作為測點(diǎn)數(shù)據(jù)，選取包含載荷作用點(diǎn)在內(nèi)的50個節(jié)點(diǎn)的位移模態(tài)信息進(jìn)行識別。

矩陣ψ的條件數(shù)為204.635 9遠(yuǎn)大于1，矩陣病態(tài)需要引入正則化參數(shù)λ。

從圖4可以看出所述的應(yīng)變位移轉(zhuǎn)化關(guān)系能夠在極小的誤差范圍內(nèi)將應(yīng)變響應(yīng)轉(zhuǎn)化為位移響應(yīng)。

圖4 無噪聲正則化時(shí)轉(zhuǎn)化過程絕對誤差Fig.4 Absolute error of conversion without noise and regularization

從圖5中可以看出當(dāng)測量的響應(yīng)信號中存在噪聲干擾時(shí)，由于矩陣的病態(tài)會導(dǎo)致結(jié)果誤差較大，影響了識別結(jié)果精度。所以需要引入正則化對噪聲進(jìn)行抑制，采用n=3的擴(kuò)展Tikhonov正則化再根據(jù)廣義交叉驗(yàn)證準(zhǔn)則確定正則化參數(shù)為0.334 8，結(jié)果如圖5所示。

對轉(zhuǎn)換關(guān)系式進(jìn)行正則化后，圖6顯示誤差穩(wěn)定在了一個很小的范圍內(nèi)，沒有出現(xiàn)失真現(xiàn)象，即通過應(yīng)變響應(yīng)的識別結(jié)果與直接通過位移響應(yīng)的識別結(jié)果基本一致，識別精度都很高。

圖5 信噪比20 dB未則化時(shí)轉(zhuǎn)化絕對誤差Fig.5 Relative error of conversion with 20 dB signal-noise ratio without regularization

圖6 信噪比20 dB正則化后轉(zhuǎn)化絕對誤差Fig.6 Relative error of conversion with 20 dB signal-noise ratio with regularization

圖7 無噪聲時(shí)基于應(yīng)變響應(yīng)的識別結(jié)果Fig.7 Noise-free identification result based on strain response

從圖8可以看出識別結(jié)果嚴(yán)重失真主要原因有兩個：

(1)速度采用差分得到，會引入高頻干擾，雖然在對應(yīng)變—位移的轉(zhuǎn)換過程進(jìn)行正則化后會降低噪聲的干擾，但是速度差分后會重新引入高頻噪聲。

(2)單點(diǎn)激勵時(shí)式(30)中FI的只有一個奇異值，本例中FI的奇異值為5.743 9×10-9,由奇異值分解理論得

(35)

由于奇異值很小導(dǎo)致噪聲誤差被放大導(dǎo)致結(jié)果失真。

圖8 信噪比20 dB時(shí)識別結(jié)果Fig.8 Result with regularization when signal-noise ratio is 20 dB

單點(diǎn)激勵時(shí)，小奇異值引起的誤差由于只有一個奇異值無法通過正則化減小，只有減小響應(yīng)中的噪聲誤差err項(xiàng)。減小響應(yīng)中的噪聲誤差只有引入濾波，本文采用八階的最小二乘平滑濾波對由應(yīng)變轉(zhuǎn)換得到的位移響應(yīng)進(jìn)行濾波后進(jìn)行識別結(jié)果，如圖9所示。

圖9 信噪比20 dB時(shí)濾波后識別結(jié)果Fig.9 Result without regularization when signal-noise ratio is 20 dB

從圖9可以看出當(dāng)對轉(zhuǎn)換得到的位移響應(yīng)進(jìn)行濾波處理后，識別結(jié)果基本處在了理論值附近，大大提高了識別進(jìn)度和識別模型的魯棒性。同時(shí)也得出當(dāng)載荷作用點(diǎn)較少時(shí)，或是系統(tǒng)矩陣雖然條件數(shù)不大但是奇異值均較小時(shí)，采用正則化的處理對噪聲的抑制作用不大，此時(shí)需要引入濾波方法對響應(yīng)進(jìn)行處理。

5 自由梁的錘擊載荷識別

為了驗(yàn)證所述方法在實(shí)際中的應(yīng)用，本文設(shè)計(jì)了如圖10所示的自由梁實(shí)驗(yàn)裝置，梁的基本參數(shù)為長1 200 mm,寬65 mm,厚4.5 mm，密度為7 850 kg/m3,彈性模量為210 000 MPa,泊松比為0.33。利用橡皮繩將梁豎直懸掛以近似自由-自由約束，在梁的一側(cè)均勻布置13個FBG傳感器記錄傳感器應(yīng)變響應(yīng)。該梁的位移模態(tài)和應(yīng)變模態(tài)信息如表1 和表2 所示。采用力錘敲擊梁未布置傳感器的一側(cè)與傳感器相對的位置，并記錄力的數(shù)據(jù)。實(shí)驗(yàn)結(jié)果如圖11所示，識別結(jié)果如圖12所示。對信號進(jìn)行濾波時(shí)仍采用平滑濾波，但是在靠近峰值時(shí)濾波效果逐漸減弱，以減輕濾波對沖擊響應(yīng)幅值的影響。

圖10 實(shí)驗(yàn)裝置圖Fig.10 Experiment equipment

特征頻率/Hz位移模態(tài)振型16.61045.77689.717148.270

表2 位移模態(tài)信息

從圖12可以看出識別出的錘擊力的大小和作用時(shí)刻與測量的力信號基本一致，作用時(shí)間由于噪聲影響較實(shí)際偏長。實(shí)際中對瞬態(tài)沖擊響應(yīng)的幅值和作用時(shí)刻更為關(guān)心，所以本文方法在實(shí)際中仍有很好的價(jià)值。

圖11 錘擊點(diǎn)應(yīng)變響應(yīng)Fig.11 Strain of the impact point

圖12 識別結(jié)果與測量得到的力信號Fig.12 Identification result and measured force

6 結(jié) 論

本文所提方法主要是基于結(jié)構(gòu)的應(yīng)變模態(tài)參數(shù)和位移模態(tài)參數(shù)，運(yùn)用應(yīng)變響應(yīng)來識別載荷。一方面由于光纖光柵的廣泛應(yīng)用，運(yùn)行飛行器的多通道動應(yīng)變響應(yīng)測量進(jìn)入實(shí)用化，因此利用多應(yīng)變完成飛行器結(jié)構(gòu)載荷識別成為可能；另一方面，通常應(yīng)變響應(yīng)對低頻載荷更為敏感，基于應(yīng)變響應(yīng)的載荷識別無疑豐富了載荷識別的手段。

本文通過仿真與實(shí)驗(yàn)算例對所提方法進(jìn)行了驗(yàn)證，結(jié)果表明：第一，方法利用應(yīng)變響應(yīng)，能夠穩(wěn)定高精度地識別結(jié)果瞬態(tài)載荷；第二，應(yīng)用位移模態(tài)和應(yīng)變模態(tài)實(shí)現(xiàn)應(yīng)變響應(yīng)到位移響應(yīng)的轉(zhuǎn)化，可以在很高的精度和穩(wěn)定性下實(shí)現(xiàn)轉(zhuǎn)換，而且過程簡單、計(jì)算量?。坏谌?，轉(zhuǎn)化過程中存在矩陣病態(tài)問題所以需要引入正則化對其進(jìn)行改善，正則化后在有噪聲干擾時(shí)能夠提高轉(zhuǎn)化精度和數(shù)值穩(wěn)定性；第四，識別過程由于采用差分求解速度會引入高頻干擾，當(dāng)系統(tǒng)矩陣維度減小時(shí)正則化對噪聲的抑制作用減小，維度為1時(shí)正則化無作用，必須采用濾波方法提高識別的精度和魯棒性。