基于快速Fourier變換法的廣義特征值問題重根辨識方法?

曹艷君 王 皓

（復(fù)旦大學(xué)航空航天系，上海 200433）

基于快速Fourier變換法的廣義特征值問題重根辨識方法?

曹艷君 王 皓?

（復(fù)旦大學(xué)航空航天系，上海 200433）

對具有重根的廣義特征值問題，采用基于快速Fourier變換的方法進行求解，實現(xiàn)重根辨識.文章中采用多次單點初始激勵的方式，仿真計算測點上的自由振動響應(yīng)，對響應(yīng)進行快速Fourier變換后得到頻域數(shù)據(jù).而后對頻域數(shù)據(jù)分析，得到固有頻率和多組測點振型數(shù)據(jù).根據(jù)單頻和重頻處的振型特性，引入振型的余弦相似度為判別參數(shù)，辨識重根.數(shù)值算例表明，該方法可有效實現(xiàn)重根辨識，同時特征值的計算能達到較高精度.

廣義特征值問題， 重根辨識， 快速Fourier變換法， 固有頻率， 動力學(xué)響應(yīng)

引言

廣義特征值問題的求解是結(jié)構(gòu)動力分析和穩(wěn)定性分析的關(guān)鍵.求解廣義特征值問題的傳統(tǒng)算法有Lanczos 法、子空間迭代法、Ritz 向量法［1，2］和 Jacobi算法等.近年來也涌現(xiàn)出一些新興算法，如混合人工魚群算法［3］、神經(jīng)網(wǎng)絡(luò)算法［4，5］、瀑布型多重網(wǎng)格法［6］等.

2014年，吳鋒、徐小明和鐘萬勰提出將快速Fourier變換法用于求解廣義特征值問題［7］.該方法的主要思想是，采用計算機仿真計算結(jié)構(gòu)的自由振動響應(yīng)，對響應(yīng)進行快速Fourier變換，從頻域數(shù)據(jù)可直觀地得到特征值信息.但對于具有重頻的廣義特征值問題，作者并未敘及.

重頻是多自由度系統(tǒng)中固有頻率重疊的現(xiàn)象，在動力學(xué)分析中較為常見.在求解廣義特征值問題過程中若沒有考慮重頻存在的可能性，求解時會遺漏重頻處的振型.這會導(dǎo)致分析結(jié)果誤差增大，甚至得到完全錯誤的結(jié)果［8］.因此，在廣義特征值問題的求解中，有必要對重頻進行辨識.

本文主要對具有重頻的廣義特征值問題進行研究，采用多次單點初始激勵，仿真計算測點上的動力學(xué)響應(yīng).在對頻域數(shù)據(jù)進行分析時，除提取固有頻率信息以外，還同時得到其對應(yīng)的振型.得到固有頻率和多組振型數(shù)據(jù)后，根據(jù)單頻和重頻處振型的不同特性，以振型向量的余弦相似度為判別參數(shù)辨識重頻.

1 重頻結(jié)構(gòu)與仿真計算

1.1 重頻結(jié)構(gòu)的自由振動問題

無阻尼系統(tǒng)的自由振動方程為:

其中x為結(jié)構(gòu)的位移向量，M為質(zhì)量矩陣，K為剛度矩陣，n為系統(tǒng)的總自由度數(shù)．

本文僅討論K為正定矩陣的情形，此時ωi大于零．

若系統(tǒng)固有頻率中存在重頻，假設(shè)ωm為k重頻（k＞1），即:

對應(yīng)的 k 個主振型記為 φm，φm+1，…，φm+k-1．系統(tǒng)的自由振動的解寫成以下形式:x（t） ＝

對于非零位移、零速度初始條件，系統(tǒng)自由振動的響應(yīng)式為:

對于第j個自由度，有響應(yīng):

式（6）中，單頻 ωs前系數(shù)為而重頻ωm前系數(shù)為在仿真計算中，這些系數(shù)可通過模態(tài)參數(shù)辨識方法得到．

第二次計算 ωs處的振型為＝，而ωm處的振型為:

由此可知，對于具有二重頻的問題，單頻處兩次計算得到的振型歸一化后是一致的．而重頻處，只要:

得到的兩個振型就線性獨立．于是，根據(jù)此特性可以對重頻進行判別．

自由振動響應(yīng)算法可采用精細積分法［9］、Runge-Kutta法［10］等．實際計算中，無需在所有的自由度上得到完整的振型．只需在結(jié)構(gòu)中選取若干個測點，得到這些測點上的振型，便可識別出重頻．

1.2 快速Fourier變換及振型提取

在利用快速Fourier變換進行仿真計算時，要求數(shù)據(jù)在時域和頻率均是離散的，且為有限長［11］．若采樣間隔為Δt，采樣點數(shù)為Nt，采樣得到的響應(yīng)數(shù)據(jù)序列為 xj（r）（r＝ 0，1，2，…，Nt-1），xj（r）做快速Fourier變換后得到的頻域數(shù)據(jù)為Xj（l），有:

其中 Aj（l）為 Xj（l）的幅值，φj（l）為相位．

對Xj（l）提取峰值便可得到結(jié)構(gòu)的固有頻率．同時，分析頻域數(shù)據(jù)可得到振型．此處的振型指在選取的測點上的振型，而非完整振型．

1.3 采樣時間間隔Δt和采樣點數(shù)Nt

仿真計算時，采樣時間間隔Δt和采樣點數(shù)Nt的選取需要滿足一定要求．

選取采樣間隔Δt之前，要對結(jié)構(gòu)的最大頻率作估計．由采樣定理可知，采樣頻率應(yīng)至少是最大頻率的兩倍．基于廣義蓋爾圓定理［12］，有最大特征值:

仿真計算時同樣需要考慮頻率分辨率．頻率分辨率指算法能將信號中兩個相近的譜峰保持分開的能力．若信號由兩個相近頻率組成，分別為ω0和ω0+Δω0．當 Δω0很小時，可能導(dǎo)致兩個峰合在一起．此時需要增加采樣點數(shù)Nt，使得頻域數(shù)據(jù)兩峰之間出現(xiàn)明顯的波谷．根據(jù)文獻［7］，采樣點數(shù)Nt可通過下式確定:

1.4 振型向量相似度度量

兩次計算在不同的點上給定初始條件，得到兩組頻率和測點上振型的數(shù)據(jù)．滿足一定條件時，可望重頻處得到的兩個振型向量線性無關(guān)．此處引入振型向量的余弦相似度刻畫振型的相似度．若有兩個 n 維向量 μ＝｛μ1，μ2，…，μn｝T和 η＝｛η1，η2，…，ηn｝T，則有μ和η夾角的余弦值為:

cosθ的值越接近于1，兩向量的相似度越高．| cosθ＝1時， μ 和 η 可以互相線性表出；cosθ＝0時，μ和η相互正交．

理論上，單頻處兩個振型向量的夾角余弦值|cosθ|＝1，而重頻處，夾角余弦值|cosθ|＜1．由此，以振型向量的余弦相似度為判別參數(shù)，可實現(xiàn)重頻的判別．

2 算例驗證與分析

2.1 算法

以具有二重根的廣義特征值問題為例，具體的算法如下:

（1）選定s個測點，確定識別指標h和頻率分辨率Δω0的值，由式（13）和式（14）確定時間間隔Δt和采樣點數(shù)Nt．

（2）取一個初始激勵點，給定初始條件．

（3）計算給定初始條件下測點的位移響應(yīng)，對響應(yīng)數(shù)據(jù)進行快速Fourier變換，得到頻域數(shù)據(jù)．

（4）從頻域數(shù)據(jù)中提取固有頻率和振型信息，振型做歸一化處理．

（5）另取一個初始激勵點，給定初始條件，執(zhí)行步驟（3）和（4），得到另一組頻率和振型數(shù)據(jù)．

（6）計算每一階頻率兩振型向量的夾角余弦值，判別重頻．

2.2 算例驗證

如圖1所示，有一四邊固支的對稱鋼板結(jié)構(gòu)，尺寸為 1m×1m×0.002m，彈性模量 E＝2×1011Pa，泊松比 μ＝0．3，密度 ρ＝ 7800kg／m3．采用矩形板單元進行網(wǎng)格劃分（板單元的每個節(jié)點上有3個自由度，分別為撓度和x方向、y方向的轉(zhuǎn)角）．如圖，劃分網(wǎng)格后共有8×8＝64個單元，81個節(jié)點．節(jié)點編號方式為:左下角節(jié)點為1號節(jié)點，沿左側(cè)邊向上編號，至左上角點為9號節(jié)點，而后以底邊左面第二個節(jié)點為10號節(jié)點，依次向上編號，以此類推，直至編至右上角點為81號節(jié)點．

圖1 鋼板結(jié)構(gòu)Fig.1 Steel plate

選擇5個節(jié)點計算響應(yīng)，節(jié)點編號為20、29、51、62和71號．第一次計算時，在22號節(jié)點撓度方向給定大小為0．05m初位移，第二次在65號節(jié)點撓度方向給定同樣大小的初位移．兩次計算初速度均為零．

取 Δω0＝0．5，h＝5，相應(yīng)的采樣時間間隔 Δt＝2×10-4s，計算響應(yīng)的數(shù)值方法為四階Runge-Kutta法．

提取頻域數(shù)據(jù)的峰值得到固有頻率，綜合多個測點的數(shù)據(jù)，得到如表1第二列和第四列所示的固有頻率結(jié)果．

根據(jù)式（10）和式（11）從頻域數(shù)據(jù)中提取振型，對于每個測點，其撓度方向計入振型．

其振型向量夾角余弦值為 cosθ1＝1．000．

夾角余弦值cosθ2＝-0．235．由此，可以判定處存在重頻．

由于篇幅原因，其余各階頻率處的振型不再一一列出，振型向量的夾角余弦值如表1第三列和第六列所示．

表1 固有頻率和振型夾角余弦值Table 1 Natural frequency and cosine values

表2 計算結(jié)果與精確解Table 2 Calculated and exact results

由表 1 可知，頻率 35.138Hz、76.550Hz、101.663Hz、133.225Hz和160.050Hz是重頻.用Matlab的eig函數(shù)計算精確解，將計算結(jié)果和精確解進行對比.表2中列出了前20階固有頻率的結(jié)果（帶?的為重頻）.由此可見，本方法能較好地實現(xiàn)重頻辨識，同時固有頻率計算達到較高的精度.

2.3 分析與討論

對于具有二重頻的問題，以上算法需要進行兩次響應(yīng)計算．理論上，只計算一次響應(yīng)也可實現(xiàn)重頻的辨識．計算響應(yīng)的初始條件為:在一個點上給定初位移，同時另一個點上給定初速度．此時重頻處實部分量得到的振型?φreal和虛部分量得到的振型?φimag線性無關(guān)．理論上可據(jù)此實現(xiàn)重頻的判別，但實際操作時，重頻辨識結(jié)果受初始條件的影響較大，效果不理想．因此，對于具有二重頻的廣義特征值問題，本文仍選擇求兩次響應(yīng)的方式．

若研究具有k重頻（k≥2）的問題，首先，可按照二重頻的方法，計算固有頻率并辨識重頻．在每階重頻處，可得兩個線性無關(guān)的振型．而后繼續(xù)選用新的初始激勵點，計算動力學(xué)響應(yīng)，并得到新的一組振型數(shù)據(jù)．若重頻處的振型不可用之前的振型線性表出，則表示該階頻率處又找到一個新的振型．繼續(xù)選用新的初始激勵點，計算響應(yīng)，直到重頻處的振型均可用之前的振型線性表出．此時，重頻處得到的線性無關(guān)的振型的數(shù)目即是該階頻率的重數(shù)．

3 結(jié)論

本文采用多次單點初始激勵的方式，從響應(yīng)的頻域數(shù)據(jù)中得到多組測點上的振型，引入振型的余弦相似度為判別參數(shù)辨識重根.以具有二重頻的問題為例進行特征值計算和重頻辨識.計算結(jié)果表明，振型向量的余弦值在單頻和重頻處有明顯差異，能夠根據(jù)余弦值有效辨識重頻，同時特征值計算可達到較高精度.

1 Li R C，Zhang L H.Convergence of the block Lanczos method for eigenvalue clusters.Numerische Mathematik，2015，131（1）:83～113

2 宮玉才，周洪偉，陳璞等.快速子空間迭代法、迭代Ritz向量法與迭代Lanczos法的比較.振動工程學(xué)報，2005，18（2）:227～232 （Gong Y C， Zhou H W， Chen P， et al.Comparative study of fast subspace iteration method，Ritz method and Lanczos method.Journal of Vibration Engineering， 2005，18（2）:227～232 （in Chinese））

3 黃華娟，周永權(quán)，韋杏瓊等.求解矩陣特征值的混合人工魚群算法.計算機工程與應(yīng)用，2010，46（6）:56～59（Huang H J， Zhou Y Q， Wei X Q， et al.Artificial fish school method for eigenvalue problems.Computer Engineering and Applications， 2010，46（6）:56～59 （in Chinese））

4 Nandy S，Sharma R，Bhattacharyya S P.Solving symmetric eigenvalue problem via genetic algorithms:Serial versus parallel implementation.Applied Soft Computing， 2011，11（5）:3946～3961

5 Tang Y，Li J P.Another neural network based approach for computing eigenvalues and eigenvectors of real skew-symmetric matrices.Computers＆Mathematics with Applications， 2010，60（5）:1385～1392

6 Chen H T，He Y H，Li Y，et al.A multigrid method for eigenvalue problems based on shifted-inverse power technique.European Journal of Mathematics， 2015，1（1）:207～228

7 吳鋒，徐小明，鐘萬勰.廣義特征值問題的快速傅里葉變換法.振動與沖擊， 2014，33（22）:67～71 （Wu F， Xu X M，Zhong W X.Fast Fourier transform method for generalized eigenvalue problems.Journal of Vibration and Shock，2014，33（22）:67～71 （in Chinese））

8 劉福林，閆維明，何浩祥.重頻結(jié)構(gòu)的振型識別及處理方法.工業(yè)建筑， 2011，41（S1）:368～372 （Liu F L， Yan W M，He H X.Mode shape identification and disposal in the structure of overlap frequencies.Industrial Construction， 2011，41（S1）:368～372 （in Chinese））

9 張繼鋒，鄧子辰，徐方暖等.一種新的改進精細直接積分法.動力學(xué)與控制學(xué)報， 2015，13（4）:241～245 （Zhang J F，Deng Z C，Xu F N，et al.A new improved precise direct integration method.Journal of Dynamics and Control，2015，13（4）:241～245 （in Chinese））

10 吳志橋，高普云，任鈞國.Runge-Kutta方法求解結(jié)構(gòu)動力學(xué)方程.系統(tǒng)仿真學(xué)報，2010，22（9）:2085～2090（Wu Z Q，Gao P Y，Ren J G.Runge-Kutta methods for time integration in structural dynamics.Journal of System Simulation， 2010，22（9）:2085～2090 （in Chinese））

11 陳恒亮，蔣勇.基于DSP的實數(shù)FFT算法研究與實現(xiàn).動力學(xué)與控制學(xué)報， 2005，3（2）:52～55 （Chen H L，Jiang Y.Design and realization of real FFT based on DSP.Journal of Dynamics and Control， 2005，3（2）:52～55 （in Chinese））

12 Nakatsukasa Y.Gerschgorin′s theorem for generalized ei-genvalue problems in the Euclidean metric.Mathematics of Computation， 2011，80（276）:2127～2142

A METHOD OF MULTIPLE-FREQUENCY IDENTIFICATION FOR GENERALIZED EIGENVALUE PROBLEMS BASED ON FAST FOURIER TRANSFORM

Cao Yanjun Wang Hao?
（Department of Aeronautics and Astronautics，F(xiàn)udan University，Shanghai 200433，China）

A method based on the fast Fourier transform was proposed to solve the generalized eigenvalue problems with multiple roots.This paper studied the dynamic responses of the measure nodes with nonzero initial condition on one point.Both natural frequencies and mode shapes were extracted from the data in the frequency domain.Responses under different initial conditions were calculated to get several sets of mode shapes.Taking the cosine similarity of mode shapes as the discriminant parameter，multiple roots were then identified.The numerical example shows that this method can identify the multiple roots efficiently，and the result reaches a high accuracy.

generalized eigenvalue problem， multi-frequency identification， fast Fourier transform， natural frequency， dynamic response

13 September 2016，revised 21 December 2016.

10.6052／1672-6553-2017-015

2016-09-13收到第1稿，2016-12-21收到修改稿.

?國家自然科學(xué)基金資助項目（11572089）

?通訊作者 E-mail:wanghao＠fudan.edu.cn

?The project supported by the National Natural Science Foundation of China（11572089）

? Corresponding author E-mail:wanghao＠fudan.edu.cn