譚博, 侯玉, 鄭華, 裴承鳴
(1.西北工業(yè)大學(xué) 數(shù)據(jù)處理中心, 陜西 西安 710072; 2.中航工業(yè)西安飛機(jī)集團(tuán)有限責(zé)任公司, 陜西 西安 710089)
結(jié)構(gòu)模態(tài)參數(shù)是結(jié)構(gòu)動力學(xué)的主要參數(shù)之一,具有簡明、直觀、物理概念清晰等優(yōu)點(diǎn)。通過系統(tǒng)的結(jié)構(gòu)模態(tài)參數(shù)可以建立相應(yīng)的數(shù)學(xué)模型,進(jìn)而根據(jù)載荷下的響應(yīng)對設(shè)計(jì)進(jìn)行必要的調(diào)整和修改。在工程應(yīng)用中,實(shí)際觀測到的信號一般是隨機(jī)過程,且具有信噪比低、模態(tài)密集、有效樣本短等特征,這樣一來,對基于觀測的模態(tài)參數(shù)估計(jì)方法就提出了近乎苛刻的要求。
矩陣束方法是20世紀(jì)80年代末、90年代初由Hua和Sarkar提出的一種模態(tài)參數(shù)估計(jì)方法,已經(jīng)在飛行器、汽車、橋梁以及建筑物的設(shè)計(jì)中有了廣泛的應(yīng)用。但是,在實(shí)際應(yīng)用當(dāng)中發(fā)現(xiàn),該方法的參數(shù)估計(jì)質(zhì)量經(jīng)常會受到諸如樣本長度、信噪比等的影響。為此,本文結(jié)合隨機(jī)減量技術(shù)對矩陣束方法進(jìn)行了改進(jìn),并應(yīng)用蒙特卡羅方法對改進(jìn)后方法的數(shù)值性能進(jìn)行了分析研究。
單模態(tài)自由度線性系統(tǒng)運(yùn)動微分方程為:
(1)
其位移響應(yīng)x(t)為:
(2)
選取適當(dāng)?shù)腁,選取響應(yīng)中x(ti)=A的點(diǎn),其中交點(diǎn)對應(yīng)的時刻為ti(i=1,2,…)。由于是線性系統(tǒng),對于自ti時刻起的響應(yīng)x(t-ti)可以看成是由ti時刻的初位移、初速度和隨機(jī)激勵引起的強(qiáng)迫響應(yīng)的三者線性疊加??梢杂?3)式表達(dá):
(3)
E[X(t)]=E[(x(0)|x(0)=A)]D(t)=AD(t)
(4)
上式表明子響應(yīng)過程X(t)可以看做一個初位移為A,初速度為零的自由振動響應(yīng)。對于多模態(tài)振動系統(tǒng),若其系統(tǒng)是線性非時變的,可看做是2個或多個單模態(tài)自由度系統(tǒng)的線性疊加。
所得子響應(yīng)過程X(t)可以表示為如下M個模態(tài)的指數(shù)函數(shù)的線性組合:
(5)
式中:s(t)為系統(tǒng)響應(yīng),n(t)為殘留的系統(tǒng)噪聲,t∈[0,T],T為最大觀測時間。對第i個模態(tài),算子si=-αi+jωi可用于表示模態(tài)的頻率和阻尼比系數(shù),其中ωi為角頻率,αi為阻尼比系數(shù)。
(5)式的離散時間形式為:
(6)
式中:zi=exp(siTs)為系統(tǒng)響應(yīng)的極點(diǎn),Ts為采樣周期,k=0,1,2,…N,N為最大采樣點(diǎn)數(shù)。
由X(t)可以構(gòu)造如下的Hankel矩陣:
(7)
式中:L為矩陣束參數(shù),通常取值位于N/4~N/3之間。
對矩陣Y進(jìn)行奇異值分解,得到:
X=UDVT
(8)
式中:U為(N-L)×(N-L)的正交矩陣;V為(L+1)×(L+1)的正交矩陣;D為(N-L)×(L+1)的對角陣,主對角線上的元素為Y的奇異值。
在已知模態(tài)個數(shù)M的情況下由D的前M個較大的非零奇異值形成如下新的矩陣D′:
(9)
D′的前M行為一個由σ1…σM構(gòu)成的M×M維的對角陣,其余行均為0。
(10)
(11)
(12)
工程中采集得到的信號總是存在噪聲影響,信噪比的大小會直接影響矩陣束方法的計(jì)算精度。本節(jié)將通過仿真信號對矩陣束方法在噪聲環(huán)境下的數(shù)值性能進(jìn)行研究。設(shè)置一個雙模態(tài)自由度系統(tǒng),模態(tài)參數(shù)分別設(shè)置為:f1=11 Hz,f2=17 Hz,d1=0.05,d2=0.03。使用矩陣束方法對加噪后的信號進(jìn)行模態(tài)參數(shù)估計(jì)運(yùn)算,并對計(jì)算結(jié)果進(jìn)行統(tǒng)計(jì),表1顯示的是1 000次計(jì)算的最終統(tǒng)計(jì)結(jié)果。
表1 不同信噪比下矩陣束方法的計(jì)算結(jié)果
由表1可以明顯看出,矩陣束方法估計(jì)所得頻率受噪聲影響較小,隨著信噪比下降沒有發(fā)生明顯變化,而估計(jì)所得阻尼比系數(shù)隨著信噪比下降有較大幅度變化,在噪聲能量和加噪前的原始信號能量相等的情況下,估計(jì)值的誤差甚至達(dá)到了真實(shí)值的2~3倍,這會影響矩陣束方法的工程應(yīng)用。
通過仿真實(shí)驗(yàn)結(jié)果可以看出,在其余條件均保持不變,僅改變信噪比的情況下,矩陣束方法所得的阻尼估計(jì)結(jié)果會產(chǎn)生較大的波動。
為減輕噪聲對矩陣束方法阻尼比估計(jì)結(jié)果的不良影響,下面將使用引進(jìn)的矩陣束方法對含有噪聲的信號進(jìn)行處理。同時加入常用的帶通濾波器降噪方法的計(jì)算結(jié)果以便于和本文方法進(jìn)行橫向比較,其統(tǒng)計(jì)結(jié)果分別如表2和表3所示。
表2 通過帶通濾波器處理后的計(jì)算結(jié)果
表3 改進(jìn)的矩陣束方法的計(jì)算結(jié)果
由表1和表2的對比可以看出,通過帶通濾波器以后,信號中噪聲的能量受到抑制,盡管計(jì)算所得的阻尼比隨著噪聲能量的增大依然有所增加,但是增加幅度明顯比未經(jīng)處理的信號得到的結(jié)果要小,在噪聲能量和加噪前原始信號能量相等的情況下,計(jì)算結(jié)果的誤差約為真實(shí)值的一半,計(jì)算精度已經(jīng)有所改善。但可以注意到,盡管阻尼比的計(jì)算精度有了明顯提高,但是計(jì)算所得的頻率同真值相比的差距反而增大,因此僅通過帶通濾波器不能很好地改善矩陣束方法的數(shù)值性能。
下面將使用本文方法對數(shù)據(jù)進(jìn)行處理,并對計(jì)算結(jié)果進(jìn)行統(tǒng)計(jì)。
由表3可以看出,本文方法的計(jì)算精度與傳統(tǒng)的矩陣束方法相比有大幅的提高,其計(jì)算頻率的最大誤差不到2%;當(dāng)原始信號能量與噪聲能量相等的時候,其計(jì)算阻尼的誤差約為20%,與帶通濾波器處理后的計(jì)算結(jié)果相比,精度有進(jìn)一步的提高。
為了進(jìn)一步研究3種方法的數(shù)值性能,下面將以阻尼比為0 dB,即噪聲和原始信號能量相等的情況為例,研究3種方法在此信噪比下1 000次計(jì)算結(jié)果的分布情況。
圖1和圖2展示的為1 000次計(jì)算結(jié)果的統(tǒng)計(jì)分布,真實(shí)頻率為11 Hz,真實(shí)阻尼比系數(shù)為0.05。圖中,梗狀圖為計(jì)算結(jié)果(阻尼或頻率)分布,實(shí)線為高斯擬合曲線,橫軸中點(diǎn)均為真實(shí)值;由圖a)至圖c)依次對應(yīng)矩陣方法、帶有濾波器的矩陣束方法以及本文方法的計(jì)算結(jié)果。
圖1 計(jì)算阻尼分布
由圖1中3種方法結(jié)果的對比可以看出,矩陣束法的計(jì)算結(jié)果整體偏離真實(shí)值較遠(yuǎn),而且分布較散,這表明每次的計(jì)算值會在一個較大的范圍內(nèi)波動,不利于通過計(jì)算確定模態(tài)的參數(shù)。而引入通濾波器后信號的計(jì)算結(jié)果分布集中,且計(jì)算的阻尼比系數(shù)與真值的差明顯縮小,矩陣束方法的數(shù)值性能得到了改善。本文方法的結(jié)算結(jié)果盡管不如第二種方法所得結(jié)果理想,但是其也較為集中分布在真實(shí)值附近,相比傳統(tǒng)的矩陣束方法,其具有更好的數(shù)值性能。
圖2 計(jì)算頻率分布
由圖2可以看出,在頻率的估計(jì)方面。矩陣束方法和引入濾波器的矩陣束方法的計(jì)算結(jié)果均與真實(shí)值有一定的差異,相比之下,帶有濾波器矩陣束方法的結(jié)果分布更加集中,略微地改善了矩陣束方法的估計(jì)頻率結(jié)果精度,但是仍然和真實(shí)值相差較大,而本文方法所得結(jié)果更為理想,其頻率的精度較前兩種方法的結(jié)果有明顯提高。
綜上所述,本文方法可以有效的提升矩陣束方法的參數(shù)估計(jì)結(jié)果的精度,與帶通濾波器這種傳統(tǒng)的降噪方法相比也具有一定的優(yōu)勢。
借助隨機(jī)減量技術(shù),本文提出了一種改進(jìn)的矩陣束方法,經(jīng)計(jì)算機(jī)仿真及蒙特卡羅統(tǒng)計(jì)分析,驗(yàn)證了本文方法在模態(tài)參數(shù)估計(jì)方面的合理性和有效性。通過提高矩陣束方法的計(jì)算結(jié)果精度,還可進(jìn)一步用于解決模態(tài)參數(shù)估計(jì)中的一些其他常見問題,如密集模態(tài)等??梢灶A(yù)見,該方法對于包括結(jié)構(gòu)顫振在內(nèi)的航空航天領(lǐng)域有著重要的工程實(shí)用意義。
參考文獻(xiàn):
[1] Thomas A J, Chard J, John E, Davies A, Francis M. Defining a Bearing Replacement Strategy Using Monte Carlo Methods[J]. International Journal of Quality & Reliability Management, 2011, 28(2): 155-167
[2] Daniel Potts, Manfred Tasche. Parameter Estimation for Nonincreasing Exponential Sums by Prony-like Methods[C]∥17thConference of the international Linear Algebra Society,Braunschweig,Germany, 2013: 1024-1039
[3] Chiu Jen Ku. Random Decrement Based Method for Modal Parameter Identification of a Dynamic System Using Acceleration Responses[J]. Journal of Wind Engineering and Industrial Aerodynamics, 2007, 95: 389-410
[4] 朱瑞可,李興源,王渝紅,趙睿,胡楠,李寬. 基于矩陣束算法的諧波和間諧波參數(shù)估計(jì)[J]. 華東電力, 2012, 40(3): 388-391
Zhu Ruike, Li Xingyuan, Wang Yuhong, Zhao Rui, Hu Nan, Li Kuan. Parameter Identification of Harmonics and Interharmonics Based on Matrix Pencil Algorithm[J]. East China Electric Power, 2012, 40(3): 388-391 (in Chinese)
[5] 莊銘杰. 蒙特卡洛仿真在移動通信中的應(yīng)用研究[J]. 計(jì)算機(jī)仿真, 2004, 21(11): 116-120
Zhuang Mingjie. Application of Monte-Carlo Simulation in Mobile Communications[J]. Computer Simulation, 2004, 21(11): 116-120 (in Chinese)
[6] 張敏,黃俐,李文雄,鐘學(xué)賦. 大型結(jié)構(gòu)模態(tài)參數(shù)識別研究[J]. 建筑科學(xué)與工程學(xué)報, 2013, 30(2): 49-54
Zhang Min, Huang Li, Li Wenxiong, Zhong Xuefu. Research on Modal Parameter Identification on Large-Scale Structure[J]. Journal of Architecture and Civil Engineering, 2013, 30(2): 49-54 (in Chinese)
[7] 徐利,鄒傳云,陳民,何毅. 基于矩陣束算法的極點(diǎn)提取分析[J]. 通信技術(shù), 2012, 45(6): 58-60
Xu Li, Zou Chuanyun, Chen Min, He Yi. Analysis on Pole Extraction Based on Matrix Pencil Algorithm[J]. Communications Technology, 2012, 45(6): 58-60 (in Chinese)