馬 超，華宏星

（上海交通大學(xué)機(jī)械系統(tǒng)與振動(dòng)國(guó)家重點(diǎn)實(shí)驗(yàn)室，上海 200240）

一種基于新的正則化技術(shù)的沖擊載荷識(shí)別法

馬 超，華宏星

（上海交通大學(xué)機(jī)械系統(tǒng)與振動(dòng)國(guó)家重點(diǎn)實(shí)驗(yàn)室，上海 200240）

針對(duì)沖擊載荷和響應(yīng)的卷積積分關(guān)系，將卷積積分離散化，從而把沖擊載荷識(shí)別問(wèn)題轉(zhuǎn)為對(duì)簡(jiǎn)單代數(shù)方程的求解。由于沖擊載荷識(shí)別問(wèn)題往往是不適定的，為了獲得穩(wěn)定解，常采用Tikhonov正則化技術(shù)來(lái)處理，可識(shí)別出的載荷的精度不是很高，提出采用一種新的正則化技術(shù)對(duì)該問(wèn)題進(jìn)行了處理。數(shù)值結(jié)果表明，該方法相對(duì)于傳統(tǒng)Tikhonov正則化方法具有更強(qiáng)的抗噪性和魯棒性，為載荷識(shí)別的高精度提供參考。

卷積積分；沖擊載荷識(shí)別；病態(tài)問(wèn)題；穩(wěn)定解；正則化

在碰撞和沖擊工程問(wèn)題中，準(zhǔn)確的獲得沖擊載荷對(duì)結(jié)構(gòu)的設(shè)計(jì)和診斷具有重要的意義。在早期的研究中利用牛頓第二定理根據(jù)結(jié)構(gòu)的質(zhì)量和實(shí)測(cè)的加速度來(lái)估計(jì)沖擊載荷。然而這種方法忽略了結(jié)構(gòu)的沖擊變形，只有假定整個(gè)結(jié)構(gòu)為剛體時(shí)才有效。接著研究人員又提出了采用力傳感器的辦法直接測(cè)量沖擊載荷，可是有時(shí)由于環(huán)境的復(fù)雜性及結(jié)構(gòu)本身的特點(diǎn)，力傳感器往往無(wú)法布置。因此，利用測(cè)量結(jié)構(gòu)的響應(yīng)來(lái)估計(jì)沖擊載荷的技術(shù)得以發(fā)展。

目前，許多研究人員對(duì)沖擊力識(shí)別技術(shù)進(jìn)行了研究，涌現(xiàn)出了很多方法。Chandrashekhara等［1］提出了基于神經(jīng)網(wǎng)絡(luò)技術(shù)對(duì)作用在復(fù)合板上沖擊載荷力重構(gòu)方法。穆騰飛等［2］提出了基于ARX模型對(duì)復(fù)合材料結(jié)構(gòu)的沖擊載荷力時(shí)間歷程的重建方法，取得了較好的效果。嚴(yán)剛等［3］提出了一種基于智能優(yōu)化算法的沖擊載荷識(shí)別的方法，該方法的基本思想是利用結(jié)構(gòu)在沖擊載荷作用下的響應(yīng)模型，通過(guò)最小化理論模型計(jì)算結(jié)果與實(shí)際測(cè)量信息之間的差別，使用微種群遺傳算法搜索出描述沖擊力時(shí)間歷程的參數(shù)，并取得了較高的精度。Bateman［4］提出采用計(jì)權(quán)加速度技術(shù)估計(jì)作用在結(jié)構(gòu)上的沖擊載荷力，但是不足之處是確定計(jì)權(quán)因子是非常困難的。Yen等［5－6］在已知沖擊位置條件下，建立了矩形板結(jié)構(gòu)中應(yīng)變響應(yīng)與沖擊載荷之間的卷積關(guān)系，并將卷積積分在時(shí)域離散成一組線(xiàn)性方程，接著定義一個(gè)數(shù)值計(jì)算結(jié)果與實(shí)際測(cè)量響應(yīng)之間的目標(biāo)函數(shù)，采用共軛梯度法求解沖擊力時(shí)間歷程。Jacquelin［7］提出采用時(shí)域反卷積技術(shù)來(lái)反演沖擊載荷。

本文在文獻(xiàn)［7］思想的基礎(chǔ)上，把沖擊載荷的識(shí)別問(wèn)題轉(zhuǎn)化為對(duì)簡(jiǎn)單的代數(shù)方程的求解。由于載荷的識(shí)別問(wèn)題往往是不適定的，為了獲得穩(wěn)定解，需采用正則化技術(shù)來(lái)處理。文獻(xiàn)［8］在理論上提出一種新的正則化算子，本文把該正則化算子引入到?jīng)_擊載荷的識(shí)別中，并把識(shí)別結(jié)果與傳統(tǒng)Tikhonov正則化技術(shù)識(shí)別結(jié)果進(jìn)行了對(duì)比研究，并討論了識(shí)別過(guò)程中測(cè)量噪聲和測(cè)點(diǎn)位置對(duì)識(shí)別結(jié)果的影響。

1 沖擊載荷識(shí)別模型的建立

在進(jìn)行沖擊載荷識(shí)別前，首先對(duì)所研究的結(jié)構(gòu)進(jìn)行以下假設(shè)：

（1）假定整個(gè)結(jié)構(gòu)為線(xiàn)性彈性系統(tǒng)。

（2）假定結(jié)構(gòu)受沖擊時(shí)變形足夠小以致于可忽略結(jié)構(gòu)的幾何非線(xiàn)性。

根據(jù)以上假定，結(jié)構(gòu)的響應(yīng)與沖擊載荷之間的關(guān)系可以表示為：

式中：h（t）表示沖擊激勵(lì)點(diǎn)和測(cè)點(diǎn)間的脈沖響應(yīng)函數(shù)，y（t）表示測(cè)點(diǎn)處的響應(yīng)，f（t）表示沖擊載荷，符號(hào)“*”表示卷積。

對(duì)式（1）進(jìn)行離散化處理，可表示為［9］：

Δt為采樣時(shí)間間隔。方程（2）即為本文進(jìn)行沖擊載荷識(shí)別的基本方程。

2 正則化技術(shù)

一般情況下，對(duì)于方程（2）采用傳統(tǒng)最小二乘法就可求解出沖擊載荷力，其解可表示為：

式中：H+為Moore-Penrose偽逆，υi為右奇異值向量，ui為左奇異值向量，si為奇異值。但是，當(dāng)滿(mǎn)足：①矩陣H的奇異值逐漸變?yōu)榱?；②矩陣H的條件數(shù)太大，即結(jié)構(gòu)矩陣H的最大奇異值和最小奇異值比值較大；上述條件之一滿(mǎn)足或同時(shí)滿(mǎn)足時(shí)，該問(wèn)題為不適定問(wèn)題［10］。為了尋求方程（2）的一組穩(wěn)定的近似解，本文引入濾波算子gλ（s）［8］：

式中：λ為正則化參數(shù)。將式（7）代入式（6）中，方程（2）的解轉(zhuǎn)化為：

從式（8）可發(fā)現(xiàn)，當(dāng)σ＝2時(shí)，上式即轉(zhuǎn)化為傳統(tǒng)的Tikhonov正則化解（Tik），因此，算子gλ（s）是包含傳統(tǒng)Tikhonov正則化算子的一種改進(jìn)正則化算子，文中將式（8）稱(chēng)為改進(jìn)正則化解（Imp）。在文獻(xiàn)［8］中指出隨著σ的增大，正則解的相對(duì)誤差的收斂階隨之提高，為了與Tik進(jìn)行對(duì)比，本文選取σ＝3來(lái)反演載荷。從式（8）中發(fā)現(xiàn)正則化參數(shù)λ對(duì)最終的解起著重要的作用，當(dāng)選取的正則參數(shù)較大時(shí)，不能很好的識(shí)別載荷；當(dāng)選取的較小時(shí)，載荷識(shí)別的正則化解將是不穩(wěn)定的，不能合理的逼近要識(shí)別的載荷。因此，合理的選取正則化參數(shù)是正則化求解成功與否的關(guān)鍵所在。目前較常用的方法有 L曲線(xiàn)準(zhǔn)則［11］和 GCV準(zhǔn)則［12］，然而，Hence［13］指出L曲線(xiàn)有時(shí)過(guò)于平滑，很難找到曲線(xiàn)上彎曲最厲害的那點(diǎn)對(duì)應(yīng)的λ值，因此本文采用GCV準(zhǔn)則來(lái)選取最佳正則化參數(shù)，GCV函數(shù)表達(dá)式為［14］：

式中：Hreg＝（HHT+λI）－1HT，并滿(mǎn)足Freg＝HregY。當(dāng)GCV函數(shù)取最小值時(shí)對(duì)應(yīng)的λ值即為最佳正則化參數(shù)。

3 數(shù)值研究及討論

3.1 算例1

圖1為一簡(jiǎn)支梁模型示意圖，其主要參數(shù)為：長(zhǎng)L＝1 m，寬W＝5 mm，高H＝10 mm，彈性模量E、泊松比ν和密度ρ分別為2．1×1011Pa，0．3和7 800 kg／m3。假定沖擊載荷作用于梁的中心位置，見(jiàn)圖（1），其作用形式如式（10）所示。圖中1和2表示測(cè)點(diǎn)位置，分別位于梁的L／10和L／2處。

圖1 簡(jiǎn)支梁模型示意圖Fig．1 The schematic diagram for the simple supporting beam

在利用式（8）進(jìn)行沖擊載荷識(shí)別時(shí)，前提條件是已知脈沖響應(yīng)函數(shù)。對(duì)于該結(jié)構(gòu)，其脈沖響應(yīng)函數(shù)理論解［15］表示為：

式中：Ycal為計(jì)算響應(yīng)值，Rnoise為正態(tài)隨機(jī)噪聲，α為噪聲因子，std（·）表示標(biāo)準(zhǔn)差，Y為含有噪聲的響應(yīng)值。最后為評(píng)估識(shí)別精度，定義量化指標(biāo)：

式中：fid為識(shí)別載荷，ftrue為真實(shí)載荷。

圖2 最小二乘法辨識(shí)結(jié)果圖Fig．2 The identification result by the least squaresmethod

表1 系數(shù)矩陣H的最大和最小奇異值Tab.1 The lowest and greatest singular value of thematrix H

圖3 Piacrd圖Fig．3 The Picard plot

圖4 GCV函數(shù)圖（a）和L曲線(xiàn)圖（b）Fig．4 GCV function plot（a）and L curve plot（b）

圖5 辨識(shí)結(jié)果圖Fig．5 The force identification result

為驗(yàn)證方法的可靠性，接下來(lái)調(diào)整噪聲因子α，對(duì)不同噪聲下的響應(yīng)數(shù)據(jù)進(jìn)行了載荷識(shí)別研究，表2列出了不同噪聲下的識(shí)別誤差值。為節(jié)省篇幅，本文只繪制了在α＝0．01條件下不同測(cè)點(diǎn)的識(shí)別結(jié)果（見(jiàn)圖6）。對(duì)表中的數(shù)據(jù)和圖進(jìn)行研究可發(fā)現(xiàn)：隨著噪聲干擾的加強(qiáng)，利用Tik和Imp法的識(shí)別精度都是逐漸減低，但識(shí)別誤差在可接受范圍內(nèi)，此外，Imp法的識(shí)別精度要高于Tik法的識(shí)別精度。這說(shuō)明Imp法相對(duì)與Tik法具有更好的抗干擾性和魯棒性。對(duì)測(cè)點(diǎn)1和測(cè)點(diǎn)2的辨識(shí)結(jié)果作進(jìn)一步分析可發(fā)現(xiàn)：測(cè)點(diǎn)1的識(shí)別精度要好于測(cè)點(diǎn)2的識(shí)別精度。結(jié)合簡(jiǎn)支梁的振型可知，這主要是測(cè)點(diǎn)2的位置正好位于結(jié)構(gòu)的第二階振型的駐點(diǎn)，而測(cè)點(diǎn)1卻避開(kāi)了該結(jié)構(gòu)前幾階振型的駐點(diǎn)，從而得出測(cè)點(diǎn)的布置應(yīng)避開(kāi)結(jié)構(gòu)前幾階振型的駐點(diǎn)，該結(jié)論與文獻(xiàn)［17］中的結(jié)論相吻合。

圖6 在α＝0．01噪聲情況下不同測(cè)點(diǎn)處辨識(shí)結(jié)果圖Fig．6 The identification results withα＝0．01 noise and differentmeasuring locations

表2 不同測(cè)點(diǎn)不同噪聲情況下作用力F的識(shí)別誤差Tab.2 The identification errors of Fwith differentmeasuring locations and noises

3.2 算例2

圖2為一加筋板結(jié)構(gòu)模型圖，板長(zhǎng)L＝540 mm，寬W＝350 mm，具體幾何參數(shù)見(jiàn)圖中所示。其物理參數(shù)為：彈性模量E＝2．1×1011Pa、泊松比ν＝0．3和密度ρ＝7 800 kg／m3。沖擊力和測(cè)點(diǎn)的具體位置分布見(jiàn)表3。待識(shí)別沖擊載荷的作用形式為：

圖7 加筋板Fig．7 Stiffened panel

表2 沖擊力和測(cè)點(diǎn)的無(wú)量綱坐標(biāo)值Tab.2 Non-dimensional positionsof forces and responses

與算例1類(lèi)似，在對(duì)沖擊載荷進(jìn)行識(shí)別前需獲得脈沖響應(yīng)函數(shù)h（t），但由于該結(jié)構(gòu)的復(fù)雜性，求解脈沖響應(yīng)的理論解比較困難，因此，需尋求另一種途徑來(lái)獲取脈沖響應(yīng)，通過(guò)研究方程（1）發(fā)現(xiàn)，該表達(dá)式可以表示為：

接著采用第一小節(jié)的類(lèi)似地離散法，上式可表示為：

因此在已知沖擊力作用的條件下，結(jié)合測(cè)量的作用力下的響應(yīng)值，就可利用第二小節(jié)的正則化技術(shù)辨識(shí)出脈沖響應(yīng)函數(shù)h（t）。本文假定已知沖擊力的表達(dá)形式為：

圖8 脈沖響應(yīng)圖Fig．8 Impulse responses

圖9 α＝0．01噪聲條件下的Piacrd圖Fig．9 The Picard plotwithα＝0．01 noise

表3 不同噪聲情況的識(shí)別誤差Tab.3 The identified errorswith different noises

圖10 α＝0．01噪聲條件下的辨識(shí)結(jié)果圖Fig．10 The identified results withα＝0．01 noise

4 結(jié) 論

針對(duì)沖擊載荷和響應(yīng)的卷積積分關(guān)系，本文將卷積積分離散化，從而把沖擊載荷識(shí)別問(wèn)題轉(zhuǎn)為對(duì)簡(jiǎn)單代數(shù)方程的求解。由于沖擊載荷識(shí)別問(wèn)題往往是不適定的，為了獲得穩(wěn)定解，文中采用一種改進(jìn)正則化技術(shù)來(lái)處理。通過(guò)數(shù)值結(jié)果表明：該方法相對(duì)于傳統(tǒng)Tikhonov正則化方法在沖擊載荷識(shí)別方面具有更強(qiáng)的抗噪性和魯棒性。

［1］Chandrashekhara K，Okafor A C，Jiang Y．Estimation of contact force on composite plates using impact-induced strain and neural networks［J］．Composites Part B：Engineering，1998，29（4）：363－370．

［2］穆騰飛，周麗，趙林虎．基于ARX模型的復(fù)合材料加筋板結(jié)構(gòu)沖擊載荷時(shí)程重建［J］．振動(dòng)工程學(xué)報(bào)，2010，26（1）：25－32．

MU Teng-fei，ZHOU Li，ZHAO Lin-hu．Impact loading time history reconstruction of stiffened composite panels using ARX model［J］．Journal of Vibration Engineering，2010，26（1）：25－32．

［3］嚴(yán)剛，周麗．加筋復(fù)合材料結(jié)構(gòu)的沖擊載荷識(shí)別［J］．航空學(xué)報(bào)，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．

［4］Bateman V，Garne T，McCall D．Force reconstruction for impact tests of an energy-absorbing nose［R］．Sandia National Labs．，Albuquerque，NM（USA），1990．

［5］Wu E，Yeh JC，Yen C S．Identification of impact forces at multiple locations on laminated plates［J］．AIAA Journal，1994，32（12）：2433－2439．

［6］Yen CS，Wu E．On the inverse problem of rectangular plates subjected to elastic impact．I：Method development and numerical verification［J］．Journal of Applied Mechanics，1995，62（3）：692－698．

［7］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．［8］范小平．不適定問(wèn)題的穩(wěn)定化算法設(shè)計(jì)及應(yīng)用［D］．淄博：山東理工大學(xué)2006．

［9］毛玉明．動(dòng)載荷反演問(wèn)題時(shí)域分析理論方法和實(shí)驗(yàn)研究［D］．大連：大連理工大學(xué)，2010．

［10］毛玉明，郭杏林，趙巖，等．自由－自由運(yùn)行體系動(dòng)態(tài)載荷反演問(wèn)題研究［J］．計(jì)算力學(xué)學(xué)報(bào)，2010，27（1）：35－39．

MAO Yu-ming，GAO Xing-lin，ZHAO Yan，et al．Study of dynamic force identification for free-free structural system［J］．Chinese Journal of Computational Mechanics，2010，27（1）：35－39．

［11］Hansen P C．The L-curve and its use in the numerical treatment of inverse problems［M］．IMM，Department of Mathematical Modelling，Technical University of Denmark，1999．

［12］Golub G H，Heath M，Wahba G．Generalized crossvalidation as a method for choosinga good ridge parameter［J］．Technometrics，1979，21（2）：215－223．

［13］Hansen P C，OLeary D P．The use of the L-curve in the regularization of discrete ill-posed problems［J］．SIAMJournal on Scientific Computing，1993，14（6）：1487－1503．

［14］Mao Y，Guo X，Zhao Y．A state space force identification method based on Markov parameters precise computation and regularization technique［J］．Journal of Sound and Vibration，2010，329（15）：3008－3019．

［15］Rao SS．機(jī)械振動(dòng)［M］．李欣業(yè)張明路．北京：清華大學(xué)出版社，2009．

［16］Hansen P C．The discrete picard condition for discrete illposed problems［J］．BIT Numerical Mathematics，1990，30（4）：658－672．

［17］智浩，文祥榮，繆龍秀，等．動(dòng)態(tài)載荷的頻域識(shí)別方法［J］．北方交通大學(xué)學(xué)報(bào)，2000，24（4）：5－10．

ZHIHao，WEN Xiang-rong，MIAO Long-xiu，et al．Dynamic loading identification in frequency domain［J］．Journal of Northerrn Jiaotong University，2000，24（4）：5－10．

Im pact force identification based on im proved regularization technique

MAChao，HUAHong-xing
（State Key Laboratory of Mechanical System and Vibration，Shanghai Jiao Tong University，Shanghai200240，China）

According to the convolution integral relationship between the impact force and response，the convolution equation was discreted and the integral equation was transfered to simple algebraic equation．Due to that the force identification problem is a well-known ill-posed problem，the solution is usually unstable．In order to obtain stable solution，an improved regularization technique was introduced．Numerical testsweremade to verify the proposed method，and the resultswere compared with those identified by the traditional Tikhonov regularizationmethod．It is shown that the proposed method gives better results than the Tikhonovmethod．

convolution integral；impact force identification；ill-posed problem；stable solution；regularization

TH212；TH213．3

A

10．13465／j．cnki．jvs．2015．12．028

2014－03－26 修改稿收到日期：2014－05－15

馬超男，博士生，1982年10月生