譜矩陣試驗(yàn)條件制定方法

吳家駒 蘇華昌 于亮

譜矩陣試驗(yàn)條件制定方法

吳家駒 蘇華昌 于亮

（北京強(qiáng)度環(huán)境研究所，北京 100076）

隨著多振動臺激勵(lì)系統(tǒng)和多自由度控制系統(tǒng)設(shè)備能力的發(fā)展，六自由度隨機(jī)振動試驗(yàn)方法已被納入美軍標(biāo)。本文專注于研究多自由度試驗(yàn)系統(tǒng)輸入條件的確定，指出兩個(gè)自由度之間振動的相干度是非常重要的參數(shù)，相干為零的獨(dú)立基礎(chǔ)激勵(lì)未必反映最嚴(yán)重工況，提出了一個(gè)廣義的六自由度振動條件的制定方法。利用滑動平均將實(shí)測信號分解成平穩(wěn)載波和時(shí)變調(diào)制波，然后多個(gè)載波信號被整合成一個(gè)時(shí)間序列，通過時(shí)平均得到歸一化的六自由度參考譜矩陣。虛擬仿真試驗(yàn)結(jié)果證明了方法的保真度。

譜密度矩陣；多振動臺試驗(yàn)；厄米特矩陣；喬列斯基分解；保真度

0 引言

多維振動試驗(yàn)技術(shù)廣泛地應(yīng)用于軍用設(shè)備的研制，多振動臺試驗(yàn)方法已被納入通用的試驗(yàn)室規(guī)范。六自由度隨機(jī)振動試驗(yàn)，在我國已成功應(yīng)用于慣性測量組合的研制試驗(yàn)，為解決隔振與捷聯(lián)的矛盾和導(dǎo)航精度的射前評估提供了環(huán)境模擬手段[1]。在最近幾年，六自由度電磁振動臺的出現(xiàn)標(biāo)志著實(shí)現(xiàn)真實(shí)多自由度隨機(jī)振動試驗(yàn)?zāi)芰Φ木薮笸黄芠2]。軍用設(shè)備環(huán)境試驗(yàn)方法標(biāo)準(zhǔn)更是從真實(shí)再現(xiàn)使用環(huán)境的角度，將原來的多振動臺試驗(yàn)方法，規(guī)范為六自由度隨機(jī)振動試驗(yàn)方法[3]。

圖1說明了實(shí)施六自由度隨機(jī)振動試驗(yàn)需要三個(gè)要素：振動臺系統(tǒng)、控制系統(tǒng)和環(huán)境的規(guī)定[4]。前兩個(gè)要素主要是涉及硬件，由于采用了“冗余結(jié)構(gòu)”的概念，使得六自由度振動環(huán)境的物理試驗(yàn)成為可能。第三個(gè)要素是短板，主要在于有關(guān)外場環(huán)境的表征、測量和試驗(yàn)條件制定方法還未規(guī)范，因?yàn)榱杂啥仍囼?yàn)包含比傳統(tǒng)單自由度試驗(yàn)更多的環(huán)境信息，尤其是還涉及運(yùn)動自由度之間的關(guān)聯(lián)程度。

最近發(fā)布的美軍標(biāo)規(guī)定，用功率譜密度矩陣來表征六自由度隨機(jī)振動環(huán)境，用相干函數(shù)來表征運(yùn)動自由度間信息關(guān)聯(lián)程度，認(rèn)為低相干度會給出保守的結(jié)果[5]。標(biāo)準(zhǔn)建議只對矩陣的下三角部分作統(tǒng)計(jì)包絡(luò)處理，然后共軛復(fù)制上三角部分，形成試驗(yàn)條件。方法保證了譜矩陣的正定性，同時(shí)自譜的包絡(luò)也會獲得保守的條件?？墒菍ψ鳛閺?fù)數(shù)的互譜實(shí)行包絡(luò)操作的正當(dāng)性，尤其是獨(dú)立即保守的論斷缺乏理論上的依據(jù)。

圖1 六自由度振動試驗(yàn)文氏圖

本文擬通過數(shù)值模擬方法，解釋相干函數(shù)的物理意義。探索運(yùn)動激勵(lì)之間的相干與力激勵(lì)之間的相干和被試結(jié)構(gòu)動力學(xué)特征的關(guān)聯(lián)程度，相干為零的獨(dú)立激勵(lì)是否一定產(chǎn)生最嚴(yán)重的響應(yīng)。在得出與之不同的結(jié)論后，提出了如何從實(shí)測信號中獲得更多相干相位信息的算法，以制定出合理的試驗(yàn)條件。首先將外場試驗(yàn)測量的多路隨機(jī)時(shí)間歷程，分解成反映振動烈度的調(diào)制波和單位高斯分布平穩(wěn)載波[6]，然后用不同方法分開處理。前者按對數(shù)正態(tài)分布統(tǒng)計(jì)，形成時(shí)間窗函數(shù)。后者整合成的長程時(shí)間序列估計(jì)出歸一化譜矩陣，與時(shí)間窗函數(shù)一起成為六自由度參考譜矩陣。為了敘述方便，將標(biāo)準(zhǔn)建議的方法稱之為“相平均法”，這里介紹的方法稱之為“時(shí)平均法”。

為了驗(yàn)證方法的合理性，設(shè)計(jì)了一個(gè)車載導(dǎo)彈運(yùn)輸振動簡化模型進(jìn)行虛擬試驗(yàn)，用蒙特卡洛法模擬多工況外場試驗(yàn)和室內(nèi)試驗(yàn)，以計(jì)算應(yīng)力估計(jì)的累積損傷作為判據(jù)，對兩種方法的效果做出定性比較。

1 參考譜密度矩陣分解

運(yùn)載器在運(yùn)動過程中受外力（風(fēng)、噪聲和路面不平等）激勵(lì)，通過安裝界面?zhèn)鬟f給設(shè)備，使之產(chǎn)生空間隨機(jī)振動，界面運(yùn)動（加速度）的大小反映出外界激勵(lì)（力）的大小。多自由度振動試驗(yàn)條件通常用頻域的譜密度矩陣來表征，其對角線項(xiàng)是自譜，反映了各自由度自身振動烈度，非對角線項(xiàng)是互譜，反映各自由度之間的時(shí)差（相位）和因果（相干）關(guān)系。試驗(yàn)?zāi)M自由度數(shù)的選擇取決于對振源的了解和室內(nèi)試驗(yàn)?zāi)M能力，隨著數(shù)據(jù)處理硬件軟件的發(fā)展和先進(jìn)控制策略的應(yīng)用，商品化的振動控制儀已經(jīng)能夠?qū)崿F(xiàn)基于譜矩陣的高斯隨機(jī)振動試驗(yàn)控制。

根據(jù)運(yùn)動控制試驗(yàn)原理，模擬了安裝界面的運(yùn)動等同于模擬了界面力[7]，這樣試驗(yàn)條件制定問題就變成如何表征界面運(yùn)動。傳統(tǒng)的單自由度試驗(yàn)，基于設(shè)備的尺度和體積遠(yuǎn)小于運(yùn)載器的假設(shè)，將安裝界面的運(yùn)動簡化成一個(gè)質(zhì)點(diǎn)的運(yùn)動。而基于多振動臺激勵(lì)的三軸振動試驗(yàn)，雖然模擬外場使用環(huán)境能力明顯提升，但仍然沒有脫離質(zhì)點(diǎn)運(yùn)動的假設(shè)。經(jīng)過多年實(shí)踐摸索，業(yè)界目前更傾向于通過六自由度試驗(yàn)來模擬軍用設(shè)備的使用環(huán)境。六自由度試驗(yàn)基于平剖面假設(shè)，即安裝界面運(yùn)動可以用剛體沿三正交軸的平移和繞三正交軸的轉(zhuǎn)動來表征。這個(gè)假設(shè)顯然受頻率范圍的限制，因?yàn)樵谠囼?yàn)頻率內(nèi)安裝界面存在彈性模態(tài)時(shí)，就需要更多自由度來表征。代表多自由度振動試驗(yàn)條件的參考譜密度矩陣，通常有兩種生成方法。一種是根據(jù)外場試驗(yàn)所測量的多路時(shí)間歷程通過譜分析統(tǒng)計(jì)形成[5]。另一種是從數(shù)據(jù)庫中獲得自譜，互譜則通過研究某些臨界響應(yīng)的所謂最小驅(qū)動方法得到[8]。兩種方法給出的參考譜矩陣都必須與控制儀控制算法匹配。

現(xiàn)代控制儀多用基于周期圖的韋爾奇方法進(jìn)行譜估計(jì)，用“逆韋爾奇”方法隨機(jī)化。試驗(yàn)?zāi)M的關(guān)鍵是要保證參考譜矩陣對每一個(gè)頻率都能實(shí)現(xiàn)喬列斯基分解，即在數(shù)學(xué)上需要矩陣滿足正定條件。以二自由度為例，譜密度矩陣可寫成

對于用單次測量時(shí)域信號變換形成的譜矩陣可以滿足正定性，但通常試驗(yàn)條件不能僅憑一次測量確定，還需要通過統(tǒng)計(jì)多次測量以預(yù)測最惡劣環(huán)境。制定試驗(yàn)條件過程中的數(shù)據(jù)處理操作，包括考慮各種工況量級的統(tǒng)計(jì)包絡(luò)、預(yù)計(jì)經(jīng)歷時(shí)間長短所做的時(shí)間壓縮等會破壞矩陣的正定性。另外，為了降低界面彈性模態(tài)影響，還需要采用冗余測量方法，通過傳遞矩陣的偽逆，來獲得反映安裝界面剛體運(yùn)動的六自由度譜密度矩陣，也會使得某些頻率不能滿足正定條件。為保證參考譜密度矩陣是厄米特陣，文獻(xiàn)[5]建議采用譜密度矩陣的下三角部分或喬列斯基分解后的下三角陣進(jìn)行加權(quán)統(tǒng)計(jì)平均（相平均），然后共軛復(fù)制上三角部分，并且基于“獨(dú)立是保守的”論斷，將出現(xiàn)小相干函數(shù)值的頻帶視作零相干，以此來簡化參考譜型。該方法能維持矩陣正定性，但其結(jié)果的保守性值需進(jìn)一步研究。下面通過數(shù)值計(jì)算，來論證零相干的獨(dú)立激勵(lì)未必產(chǎn)生極值響應(yīng)，并研究如何處理相干函數(shù)，才能使參考譜矩陣正定且結(jié)果保守。

2 相干函數(shù)影響分析

或者

左端項(xiàng)代表外場測量數(shù)據(jù)估計(jì)的響應(yīng)譜矩陣，其對角線項(xiàng)包含的相干函數(shù)與系統(tǒng)的頻響特征和激勵(lì)矩陣關(guān)聯(lián)，其解析形式比較復(fù)雜，推導(dǎo)過程繁瑣。

本文擬通過數(shù)值分析方法，定性地解釋相干函數(shù)的物理意義。頻域乘積可轉(zhuǎn)化為時(shí)域的卷積，因此響應(yīng)的時(shí)間歷程可表示為

將圖2的模型進(jìn)一步簡化，用一根剛性梁兩點(diǎn)彈性支撐的二自由度系統(tǒng)表示，以便進(jìn)行兩種簡化梁的數(shù)值計(jì)算仿真。首先設(shè)定振源為獨(dú)立白噪聲，激勵(lì)譜密度矩陣為單位陣，即互譜為零，自譜為1進(jìn)行仿真。圖3給出二自由度系統(tǒng)（剛性梁）響應(yīng)相干隨頻率變化（較高位置）曲線，圖中還給出了系統(tǒng)的頻響交聯(lián)項(xiàng)（較低位置）?？梢钥闯?，在低頻相干度接近于1，僅在（跨點(diǎn)）頻響的谷值處出現(xiàn)低值（約為0.1）。圖4給出彈性梁（多自由度系統(tǒng)）的響應(yīng)相干隨頻率變化曲線?？梢钥闯觯冢琰c(diǎn)）頻響谷值頻帶內(nèi)出現(xiàn)低的相干度，具體位置與此頻帶內(nèi)的反共振谷深度有關(guān)。圖3和圖4說明響應(yīng)相干與系統(tǒng)的動力學(xué)特征有關(guān)。其次設(shè)定激勵(lì)譜密度矩陣為對角陣，但自譜量級不同進(jìn)行仿真。按相同方法分析結(jié)果見圖5?？梢钥闯觯宰V量級變化會改變出現(xiàn)最小相干函數(shù)值的頻率。最后，設(shè)定兩個(gè)激勵(lì)的自譜密度相同但相干不同進(jìn)行仿真，按相同方法分析結(jié)果見圖6?？梢钥闯?，響應(yīng)相干函數(shù)值與激勵(lì)相干函數(shù)值相應(yīng)變化趨勢相同。不過，響應(yīng)相干出現(xiàn)最低值的頻率也發(fā)生變化，而且即使是全相干激勵(lì)，在某些頻率響應(yīng)的相干函數(shù)也不為1。

通過數(shù)值仿真分析可以看到，外場測量數(shù)據(jù)本質(zhì)上是載體在外力作用下的運(yùn)動響應(yīng)。各個(gè)測量自由度之間的相干與相位，自然與外力之間的相對烈度和相干有關(guān)，也與振動系統(tǒng)固有特征有關(guān)。自譜的大小反映振動的烈度，可以通過提高自譜的量級來增加試驗(yàn)的保守性。但是相干函數(shù)能否給出一個(gè)保守性的指標(biāo)，還需要深入探討。

圖3 二自由度系統(tǒng)響應(yīng)相干函數(shù)

圖4 多自由度系統(tǒng)響應(yīng)相干函數(shù)

圖5 相干函數(shù)隨激勵(lì)量級比變化

圖6 相干函數(shù)隨激勵(lì)相干變化

進(jìn)一步分析輸入信號之間相位關(guān)系對輸出的影響，有助于外場真實(shí)環(huán)境的模擬。仍以二自由度模型為例，回到（9）式

3 譜矩陣估計(jì)的相平均法

圖7 輸入相干效果（力激勵(lì)）

圖8 輸入相干效果（運(yùn)動激勵(lì)）

圖9 相平均制定試驗(yàn)條件流程

這里的“上限包絡(luò)”指的是，多次測量得到的譜矩陣是變化的，對每一個(gè)頻率計(jì)算譜矩陣的每一個(gè)元素的均值和方差，并按一定的分布（例如對數(shù)正態(tài)分布），估計(jì)最嚴(yán)重環(huán)境，作為試驗(yàn)條件?？墒腔プV是復(fù)數(shù)，實(shí)部和虛部有正有負(fù)，不會與自譜有相同的分布。實(shí)踐證明互譜的平均結(jié)果是趨于降低相干度。雖然互譜還能用極坐標(biāo)表示（模和相位），但是由于周期性纏繞，對相位進(jìn)行求均值運(yùn)算會得出不正確的結(jié)果。通過以上分析發(fā)現(xiàn)，相平均方法的保真度可能存在問題。既然運(yùn)動激勵(lì)的相干和相位對響應(yīng)有明顯影響，那么就有必要在確定參考譜矩陣時(shí)，盡量維持原有的互譜屬性。為了實(shí)現(xiàn)這一目的，本文采用區(qū)別對待矩陣對角線項(xiàng)和非對角線項(xiàng)的制定試驗(yàn)條件方法。

4 譜矩陣估計(jì)的時(shí)平均法

首先將外場測量的多路隨機(jī)時(shí)間序列，逐一分解成一個(gè)窗函數(shù)（調(diào)制波）和單位方差的高斯分布隨機(jī)序列（載波），然后用不同的方法統(tǒng)計(jì)處理[6]。對于一次測量

5 虛擬仿真試驗(yàn)

為了研究方法的保真度，繼續(xù)利用簡化版的車載導(dǎo)彈運(yùn)輸振動模型，進(jìn)行兩點(diǎn)激勵(lì)的均勻梁虛擬試驗(yàn)，用蒙特卡洛法模擬多工況外場試驗(yàn)和室內(nèi)試驗(yàn)，以計(jì)算應(yīng)力估計(jì)的累積損傷作為判據(jù)，對不同方法的效果做出定性比較。

一共完成了九次試驗(yàn)，每一次先產(chǎn)生烈度（均方根值）隨機(jī)設(shè)置的相互獨(dú)立的白噪聲，作為梁的外激勵(lì)力。計(jì)算兩個(gè)激勵(lì)的加速度響應(yīng)時(shí)間歷程，模擬一個(gè)工況的外場測量，作為制定參考譜密度矩陣的一個(gè)子樣；同時(shí)計(jì)算出梁的應(yīng)力時(shí)間歷程，并估計(jì)所造成的累積損傷。然后，用九個(gè)工況激勵(lì)點(diǎn)的加速度時(shí)間歷程，分別按相平均法和時(shí)平均法形成參考譜密度矩陣。最后，用兩個(gè)譜密度矩陣，反演成兩路加速度時(shí)間歷程，作為基礎(chǔ)激勵(lì)，計(jì)算梁的應(yīng)力響應(yīng)時(shí)間歷程和所造成的累積損傷。室內(nèi)試驗(yàn)?zāi)M也進(jìn)行了九次。

圖10 時(shí)平均制定試驗(yàn)條件流程

圖11 虛擬試驗(yàn)流程

圖12展現(xiàn)出虛擬試驗(yàn)的相對累積損傷。最淺顏色曲線表示九次實(shí)測振動激勵(lì)造成的累積損傷的上限極值。最高位置曲線和最低位置曲線分別表示基于九次實(shí)測振動，用時(shí)平均法和相平均法得到參考譜矩陣激勵(lì)造成的累積損傷的預(yù)測值下限。

可以看到，兩種方法相比，時(shí)平均法的預(yù)測效果更接近外場環(huán)境，而且偏于保守。相平均法的預(yù)測效果偏離外場環(huán)境，特別是會造成“欠”試驗(yàn)，原因可用圖13來說明，圖13是運(yùn)動激勵(lì)的參考譜矩陣。

兩者自譜的全頻段幾乎一樣，兩者的相干和相位僅在300～500Hz帶內(nèi)一致，其余頻帶差別較大，且相位基本反相。

時(shí)平均法得到的相干很高，相平均法所得到的相干卻很低，而低相干趨于獨(dú)立的運(yùn)動激勵(lì)的結(jié)果不是保守的。

圖12 累積損傷結(jié)果比較

圖13 參考譜矩陣比較

6 結(jié)論

隨著振動臺和控制儀的進(jìn)步，六自由度隨機(jī)振動試驗(yàn)已經(jīng)標(biāo)準(zhǔn)化，一個(gè)新的參數(shù)-譜矩陣被用來表征多維隨機(jī)振動環(huán)境，試驗(yàn)條件的制定迎來了新挑戰(zhàn)。傳統(tǒng)基于實(shí)測外場測量的（自）譜密度統(tǒng)計(jì)方法，是否適用相干相位的統(tǒng)計(jì)處理，是一個(gè)需要討論的問題。

通過數(shù)值分析表明，制定試驗(yàn)條件的外場測量數(shù)據(jù)本質(zhì)上是載體在外力作用下的運(yùn)動響應(yīng)，各個(gè)測量自由度之間的相干相位，與外力之間的相對烈度和相干有關(guān)，也與振動系統(tǒng)的固有特征有關(guān)。由于相干與試驗(yàn)保守性之間沒有烈度那樣的相隨關(guān)系，依據(jù)外場測量數(shù)據(jù)制定譜矩陣規(guī)范時(shí)，應(yīng)當(dāng)盡可能保留特定工況下的相干特征。從這個(gè)意義上講，譜矩陣條件只能是“量身定做”，不能成為“貨架產(chǎn)品”。軍用規(guī)范推薦的相平均方法趨向于產(chǎn)生獨(dú)立的譜矩陣試驗(yàn)條件，但獨(dú)立運(yùn)動激勵(lì)未必產(chǎn)生保守的響應(yīng)。本文提出的實(shí)測信號分解→整合單位高斯載波→通過時(shí)平均得到歸一化參考譜矩陣的方法，有效地保留自由度之間相干特征。基于兩點(diǎn)激勵(lì)的車載導(dǎo)彈運(yùn)輸振動簡化模型仿真結(jié)果表明，該方法既模擬了外場環(huán)境又能做到適度保守。雖然該方法通過了虛擬試驗(yàn)演示驗(yàn)證，但仍需要更多的六自由度外場測量和物理試驗(yàn)來證明，尤期需要建議一個(gè)有效評價(jià)體系，相關(guān)問題會在后續(xù)的文章中討論。

[1] 吳家駒, 等.用于慣性測量組合的金屬減震器[R].航天報(bào)告, 2002(2): 1-5.

[2] Hoksbergen J.Bringing broadband 6-DOF field vibration environments into the lab-tensor 18kN vibration test system [J].Sound & Vibration, Dynamic Testing, March, 2013.

[3] MIL-STD-810H .Multi Exciter Test, Method 527[S].

[4] Maitin L, Schctral S.On the nature of spectral density matrices used to the characterize vibration environments[C].Proceedings of the 87thShock & Vibration Symposium.

[5] Hale M.Spectral density matrix transformation[J].Journal of the IEST, 2017, 60(1):17-26.

[6] 吳家駒, 付瑋, 張鵬飛.基于β分布隨機(jī)數(shù)排序的非高斯振動模擬方法[J].強(qiáng)度與環(huán)境, 2017, 44(2): 10-16.[Wu Jiaju, Fu Wei, Zhang Pengfei.Simulation method of non-Gaussian vibration base on re-arrangement of β distribution random number [J].Structure & Environment Engineering, 2017, 44(2): 10-16.]

[7] 吳家駒, 榮克林.多維振動環(huán)境試驗(yàn)方法[J].導(dǎo)彈與航天運(yùn)載技術(shù), 2003(4): 27-32.[Wu Jiaju, Rong Kelin.Multi-dimentional vibration environmental test method[J].Missiles and Space Vehicles, 2003(4): 27-32.]

[8] Smallwood D O.Extreme input/outputs for multiple input multiple outputs linear systems[S].SANDE 2005-5307.

[9] Smallwood D O.Multiple shaker random vibration control an update PIES[S].1999.

[10] 吳家駒, 丁鎮(zhèn)軍, 丁富海.非高斯多維隨機(jī)振動分析[J].強(qiáng)度與環(huán)境, 2018, 45(4): 1-6.[Wu Jiaju, Ding Zhenjun, Ding Fuhai.Analysis of non-Gaussian multi-random vibration [J].Structure & environment engineering,2018, 45(4): 1-6.]

[11] Sandia National Laboratories.Sandia’s new vibration table promises different ways of testing[S].2015.

[12] NASA-HDBK-7005.Dynamic environmental criteria [S].2001.

On the Derivation method of Spectral Density Matrix Specifications

WU Jia-ju SU Hua-chang YU liang

(Beijing Institute of Structure and Environment Engineering, Beijing 100076, China)

With the development of multiple excitation systems and MDOF vibration control systems, the 6-DOF laboratories vibration test method is brought into MIL-STD.This paper concentrates on the determination of an input specification for such MDOF systems.It is pointed out that the coherent value between two degrees of freedom is a very important parameter.The motion excitation of independent foundation with zero coherence does not necessarily reflect the most serious condition.A generalized 6-DOF vibration specification development technique is proposed.The measured signal is decomposed to the stationary carrier and time-varying modulated wave by the moving average.Then, a lot of carrier is integrated into a time series.The normalized 6-DOF spectral density matrix is obtained by time averaging.Virtual simulation results demonstrate the fidelity of this method.

Spectral Density Matrix (SDM);Multiple-Exciter Test(MET); Hermitian Matrix; Chlesky Decompose; Fidelity

V416.2

A

1006-3919(2021)06-0050-09

10.19447/j.cnki.11-1773/v.2021.06.008

2021-06-29；

2021-09-15

吳家駒（1939－），男，研究員，研究方向：隨機(jī)振動分析與試驗(yàn)；（100076）北京9200信箱72分箱.