基于MICCD的時變結構模態參數辨識

王 豪，藍 鯤，夏國江，耿勝男

(北京宇航系統工程研究所， 北京 100076)

1 引言

在航天領域，姿態控制系統的設計、載荷的計算、結構的優化以及元器件的安裝等都需要準確的結構模態數據。一般先通過地面模態試驗對理論模型進行修正，再利用修正的模型計算各個狀態下的模態特性。然而地面模態試驗大多局限于時不變結構，無法真實地模擬火箭飛行中的狀態，因此開展火箭飛行過程中的時變模態參數辨識對整個火箭的安全性與穩定性具有重要意義。

傳統的時變結構模態辨識通常采用“時間凍結”的理論，認為在凍結的時間窗口中結構參數保持不變，再利用成熟的時不變模態辨識方法對各個時間段進行處理最終得到時變模態參數，例如董嚴等近似地認為火箭結構的模態在每個時間段內保持不變，利用從火箭發射至發動機分離時間段內各測點的測試數據，基于自回歸滑動平均模型成功得到了各階模態參數隨時間的變化關系。王鵬輝等采用基于自然激勵法的組合方法對火箭模態參數進行識別，將不同時刻的參數進行擬合處理得到火箭的時變模態參數，獲得了模態參數隨加液氫注量增加的變化規律。但“時間凍結”法忽略了時間變化對結構參數的影響，模態辨識準確度并不高。

近年來隨著信號處理技術的發展，將信號分解技術和時頻變換相結合的時變模態辨識方法發展迅速。時頻變換可以提供頻率分辨率，信號分解可以提供時間精度，兩者的結合提高了模態辨識的準確程度。例如王佐才等將解析模態分解(analytical modal decomposition,AMD)與小波變換結合用于時變模態辨識，取小波脊線間的平均值作為AMD分解的時變截止頻率，提高了瞬時頻率的識別精度。王超等利用變分模態分解(variational mode decomposition,VMD)加廣義Morse小波的方法，對小車經過主梁時的加速度響應進行分析，成功識別了第一階模態的時變頻率。Wei等將固有啁啾分量分解方法(intrinsic chirp component decomposition,ICCD)與廣義參數化時頻變換(general parameterized time-frequency transform,GPTFT)結合應用到時變系統自由響應下的模態辨識中，準確度遠高于傳統的希爾伯特-黃變換(hilbert-huang transform,HHT)。

但上述方法始終是單通道的模態辨識方法，一次僅能處理一個測點的數據，而火箭上的測點往往是成百上千個，逐點處理不僅效率低下，也沒有充分利用各測量通道之間的相關性。鑒于目前鮮有學者開展該方面的研究，本文提出了基于多通道固有啁啾分量分解(multivariate ICCD,MICCD)的時變結構模態辨識方法，該方法不僅利用GPTFT來提高瞬時頻率(instantaneous frequency,IF)的分辨率，還通過擴展的最小二乘法求解多通道線性方程組來精確提取各時變分量，提高了時變模態參數辨識的準確性。

2 模態參數辨識

一個自由度的時變結構振動微分方程可以表示為：

(1)

式中：()、()和()分別為時變結構時刻的質量矩陣、阻尼矩陣以及剛度矩陣；()為結構的位移響應向量；()為外部激勵。根據模態疊加法，第自由度的脈沖位移響應()可以表示為各階模態響應的疊加，即：

(2)

其中

式中：為頻率分量總數；,()、,()分別為第自由度第階位移振型向量值和瞬時振幅(instantaneous amplitude,IA)；()為第階歸一化后的模態坐標；和代表()的振幅和初始相位；()和()代表時刻第階的無阻尼固有頻率和阻尼比；()則為第階的有阻尼固有頻率，由于結構的阻尼比通常很小，因此可近似認為()≈()。

值得注意的是，模態疊加法是基于時不變系統提出的。在時變系統領域，并沒有關于模態響應的精確的理論推導，大多數研究采取了時變系統響應由多分量瞬時模態響應疊加而成的假設，幸運的是，這些論文的結果證明這種假設對于時變系統是適用的，本文的仿真分析也證明了這一點。

若各通道的各階的瞬時振幅已知，則可計算出時變結構的第階歸一化瞬時振型向量()：

(3)

因此，瞬時頻率(IF)和瞬時振幅(IA)的準確提取是模態參數辨識的關鍵。

3 多通道固有啁啾分量分解

3.1 多通道固有啁啾分量

在語音處理、雷達應用、機械故障診斷等各種應用中，信號通常是非平穩的，可以被建模為調幅和調頻信號，也稱為啁啾信號。對于多通道的啁啾信號而言，以結構響應為例，多個通道的響應信號都包括了同樣的頻率分量，故定義多通道固有啁啾分量(multivariate intrinsic chirp component,MICC)如下：

(4)

式中，()即為多通道固有啁啾分量，是相同頻率分量成分的集合。結合式(2)，多通道的測量信號()可建模如下：

(5)

式中，,()=,()cos，,()=-,()sin；()則代表分解誤差和測量噪聲。與單通道ICCD方法類似，可利用冗余傅里葉級數對瞬時振幅,()與瞬時頻率()進行擬合：

(6)

其中的待求系數為

其中，=2π=2π，為信號樣本個數，當≥2時，上式變為冗余傅里葉級數模型，根據過完備字典下的信號稀疏分解理論，這樣做有利于減小擬合誤差，提高系數求解和建模的準確性，本文取=2；和分別為IA與IF的冗余傅里葉模型階數，決定了擬合的復雜程度，模型階數的確定方法可參考文獻[7]。總的來看，式(5)和式(6)說明了瞬時振幅,()與瞬時頻率()之間的線性關系，即在瞬時頻率()已知的情況下，通過求解一組線性方程組，可以很容易地得到瞬時振幅,()。

3.2 瞬時頻率提取

瞬時頻率的準確提取依賴于能量高度集中的時頻表示，能量越集中，時頻脊線提取的準確度也就越高。傳統的時頻表示方法包括短時傅里葉變換、小波變換和Wigner-Ville分布等。雖然短時傅里葉變換和小波變換在變換基的使用上有所不同，但它們本質上都是使用水平線來近似表達給定信號在時頻平面上的瞬時頻率，Wigner-Ville分布雖然可以提供高濃度的時頻表示，但在分析多分量信號時會存在交叉項干擾，頻率分辨性不強，而由Yang等提出的GPTFT在核函數表達式合理的情況下，能夠準確表征各種非平穩信號的時頻特征。為與式(6)對應，本文采用如下的核函數與變換框架：

(7)

其中

()=()+{()}

(8)

其中

3.3 瞬時振幅估計

在單通道ICCD算法中，在已知IF信息的情況下，可建立具有范數約束的最小二乘模型求解瞬時振幅(IA)。然而，該模型不能應用于多通道信號，因為多通道信號處理需要同步分解所有通道信號，并不能用向量范數逐個描述，因此需要將最小二乘模型推廣到多通道模型。為了更好地解釋MICCD方法，將多通道輸入信號模型()重寫為

=+

(9)

式中，=[,,…,]；=[,,…,]為包含模態信息的核函數，其中

= diag[cos(()) … cos((-1))]

= diag[sin(()) … sin((-1))]

為待求IA冗余傅里葉模型系數矩陣：

為使式(9)的分解誤差達到最小，同時避免求矩陣偽逆帶來的“病態”問題，引入正則化因子：

(10)

(11)

(12)

(13)

綜上所述，本文算法分為兩步：一是利用GPTFT估計多通道固有啁啾分量的瞬時頻率信息；二是利用式(11)計算冗余傅里葉級數模型的系數，然后分別用式(12)和式(13)重建瞬時振幅和分量。

4 仿真校驗

4.1 仿真設置

在Matlab/Simulink中對如圖1所示的三自由度時變結構進行仿真模擬，給定第2個自由度以白噪聲激勵，使整個系統做隨機振動，采樣頻率為100 Hz，采樣時間為15 s，求解器選擇Runge-Kutta算法。

圖1 三自由度彈簧阻尼系統Fig.1 Three-degree-of-freedom spring damping system

圖1中各項時變的物理參數為

()=05,()=15e-005,()=25(kg)

()=02,()=02

()=03+003sin(05π)

()=15 000,()=12 000

為衡量瞬時固有頻率辨識的精確度，定義如下的頻率辨識誤差指標：

(14)

模態保證準則(modal assurance criterion,MAC)能夠反映2個振型向量之間的相關程度，因此可以利用MAC值來表征某一時刻瞬時振型辨識的準確程度：

4.2 隨機振動分析

表1 各個頻率成分迭代擬合的頻率辨識誤差Table 1Frequency identification error of iterative fitting of each frequency component Hz

以提取第3階頻率分量為例，GPTFT迭代效果如圖2所示，圖2中的曲線即為擬合IF曲線，可以隨著迭代次數的增加，第3階的能量分布越來越集中，證明了GPTFT在面對頻率變化較快的分量時依然能夠提供能量高度集中的時頻分布。最終得到的各階瞬時頻率曲線如圖3，IF擬合值與理論值基本吻合，說明本算法可以準確地提取瞬時頻率。

將上述的IF信息代入式(11)～式(13)得各通道各分量的分解結果及殘差如圖4所示，圖4中,代表了第通道第階啁啾分量；代表第通道信號的分解誤差，相比原始信號分解誤差能量極小。

圖2 第3階頻率分量的擬合迭代過程示意圖Fig.2 Fitting iterative process of the 3rd frequency component

圖3 基于MICCD的瞬時頻率曲線Fig.3 Instantaneous frequency identification result based on MICCD

圖4 MICCD分解結果及線差圖Fig.4 MICCD decomposition results

圖4證明了MICCD方法不僅同時處理多個通道的信號，還可以實現同一種頻率分量的對齊，且每個分量的重構精度極高。將分解得到的瞬時振幅(IA)代入式(3)可以得到瞬時振型，以第2階振型為例，辨識振型與理論振型隨時間變化趨勢如圖5所示，可以在誤差允許的范圍內，辨識振型的變化趨勢與理論振型保持一致，說明本文方法可以有效地提取結構瞬時振型。各階振型的MAC曲線如圖6所示，任意時刻的MAC值都十分接近1，說明瞬時振型的辨識精度很高，與理論振型有很高的相關性。

圖5 第2階瞬時振型辨識結果圖Fig.5 2nd instantaneous mode shape identification result

圖6 基于MICCD的瞬時振型辨識結果曲線Fig.6 Instantaneous mode shape identification result based on MICCD

4.3 對比分析

為了展示本文方法的優越性，又引入MVMD方法對時變結構模態參數進行辨識。MVMD是在2019年由Rehman等基于單通道VMD方法提出的多通道信號分解方法，通過建立約束變分模型，實現多通道信號非遞歸的自適應分解，將多通道信號在相同的頻率尺度分解為相同數量的固有模態函數(intrinsic mode functions,IMFs)之和，保證了多通道信號分解時各階分量頻率的一致性。

在反復進行模態分解試驗后，選擇最優的分解參數，令分解的IMF數量為3，帶寬約束參數為700，迭代收斂容差為10。經MVMD分解后，每個自由度的響應均被分解為3個頻率一致的IMF分量，每一個IMF分量都反映了時變結構的某一階固有頻率。為獲得瞬時固有頻率，本文采用Hilbert變換求解每個IMF分量,的瞬時頻率和瞬時振幅：

(16)

式中:,、,和,分別為第通道第個IMF的瞬時相位、瞬時頻率與瞬時振幅。圖7展示了基于MVMD方法得到的瞬時頻率曲線，由于Hilbert方法直接求得的IF曲線呈現反復振蕩的特點，并不符合真實結構的頻率變化規律，利用式(6)對其擬合獲得擬合IF曲線，結果表明前兩階的瞬時頻率比較準確，但第3階出現了明顯的模態混疊現象，瞬時頻率曲線后半段跳躍到第2階。對比圖3發現基于MICCD方法的瞬時頻率辨識精度更高，瞬時振型亦是如此(對比圖6與圖8)。

圖7 基于MVMD的瞬時頻率辨識結果曲線Fig.7 Instantaneous frequency identification result based on MVMD

圖8 基于MVMD的瞬時振型辨識結果曲線Fig.8 Instantaneous mode shape identification result based on MVMD

5 結論

1) 引入廣義參數化時頻變換(GPTFT)，提高了時頻面內的能量集中程度，再迭代優化提高了瞬時頻率辨識的準確率；

2) 擴展最小二乘法求解多元線性方程組，充分利用了各通道間的關系，使得啁啾分量重構更加精確，獲得了精確的振幅包絡和瞬時振型；

3) MICCD方法需要對各階的瞬時頻率進行準確提取作為輸入信息，比MVMD自適應方法適用范圍更廣，分解誤差更小，模態辨識精度更高。