基于拉普拉斯特征映射的三維結構模態(tài)分析

符偉華，王 成，2+，陳建偉，賴雄鳴，李海波,2

(1.華僑大學 計算機科學與技術學院，福建 廈門 361021；2.廈門市企業(yè)互操作與商務智能工程技術研究中心，福建 廈門 361021；3.圣地亞哥州立大學 數(shù)學與統(tǒng)計學院，圣地亞哥 加州 92182，美國；

0 引言

模態(tài)是結構本身具有的振動特性，通過工作模態(tài)分析(Operational Modal Analysis, OMA)方法識別每階模態(tài)的參數(shù)(如模態(tài)固有頻率、振型、阻尼比等)，可以了解結構的動力學特性，進而進行結構的損傷識別[1]、結構優(yōu)化[2]等。然而，對于許多大型復雜結構，可以利用的激勵方式只有工作狀態(tài)下的環(huán)境激勵，從而導致無法對激勵輸入進行測量。為了解決該難題，一種僅通過利用結構的輸出響應識別結構模態(tài)參數(shù)的方法引起了許多研究者的關注[3]。

當前，工作模態(tài)分析已經得到了廣泛應用，各類工作模態(tài)識別方法層出不窮。PONCELET等[4]采用了獨立分量分析和二階盲辨識兩種盲源分離技術，進行了純輸出模態(tài)分析。王成等[5]利用了主成分分析(Principal Component Analysis, PCA)進行模態(tài)參數(shù)識別,解決其他方法識別過程中出現(xiàn)虛假模態(tài)的問題。DEVRIENDT等[6]提出的基于傳遞率測量的OMA技術，用于降低由于諧波存在而導致模態(tài)參數(shù)識別錯誤的風險。

這些模態(tài)識別方法都取得了很好的實驗效果，但對復雜三維連續(xù)體結構工作模態(tài)參數(shù)識別的研究甚少。白俊卿等[7]利用局部線性嵌入的流形學習方法，將模態(tài)數(shù)據(jù)集視作高維數(shù)據(jù)，通過尋找低維的線性嵌入，從而識別復雜三維連續(xù)體結構工作模態(tài)參數(shù)。王成等[8]利用模態(tài)坐標的概念，找出動力機械系統(tǒng)的振動模態(tài)振型與線性混合矩陣之間、模態(tài)響應與主成分之間的一一對應關系，提出一種新的三維結構時域統(tǒng)計模態(tài)分析方法。通過最小二乘反迭代，王建英[9-10]提出了二階盲辨識的三維連續(xù)體結構工作模態(tài)參數(shù)識別方法。通過振動響應三維矩陣的直接組裝，張?zhí)焓鎇11-12]提出了主元抽取的三維連續(xù)體結構工作模態(tài)參數(shù)識別方法。

近年來，流形學習在圖像處理等多個工程領域得到了廣泛應用。流形是一個幾何概念，它表示鑲嵌在高維空間中的低維幾何結構[13]。基于這個概念，流形學習的降維過程中需保持降維之后的數(shù)據(jù)同樣滿足與高維空間流形有關的幾何約束關系。目前，流形學習[14]運用較為廣泛的算法有等距離映射(Isometric feature mapping, Isomap)[15-16]、局部線性嵌入算法[17]等。拉普拉斯特征映射(Laplacian Eigenmaps, LE)[18]是流形學習的主要算法之一，其基本思想是從局部角度構建數(shù)據(jù)之間的聯(lián)系，使得相互間有聯(lián)系的點在完成降維后盡量保持該聯(lián)系，這使得降維后的數(shù)據(jù)集依然能保持原來數(shù)據(jù)集的數(shù)據(jù)結構。相比主成分分析，LE能考慮數(shù)據(jù)的內在幾何；相比等距離映射，拉普拉斯特征能僅從局部構建數(shù)據(jù)的內在幾何，使得算法時間復雜度更低[19]。本文的主要創(chuàng)新點如下：

(1)統(tǒng)一提出了基于流形學習的工作模態(tài)參數(shù)識別方法體系。其基本思想是找出各階模態(tài)坐標響應與低維嵌入數(shù)據(jù)之間的一一對應關系，從而將工作模態(tài)參數(shù)識別問題轉化為結構振動響應數(shù)據(jù)的流形降維問題。

(2)將LE算法引入動力學系統(tǒng)，提出基于LE的工作模態(tài)參數(shù)識別方法，并與主成分分析、等距離映射進行理論比較。

(3)結合最小二乘廣義逆反迭代或者直接組裝三維位移響應求解法，提出了基于拉普拉斯特征映射的復雜三維結構工作模態(tài)參數(shù)識別方法。

(4)在三維圓柱殼結構上的振動響應仿真數(shù)據(jù)集上進行基于LE的三維結構工作模態(tài)參數(shù)識別方法的效果驗證。

1 基于LE的復雜三維連續(xù)體結構模態(tài)參數(shù)識別

1.1 工作模態(tài)參數(shù)識別問題的描述

工程結構系統(tǒng)一般可看作是一個具有一定粘性比例阻尼的系統(tǒng)。結合結構的動力學理論，對于包含n個自由度的線性時不變振動系統(tǒng)，該系統(tǒng)在物理坐標系統(tǒng)中的動力學方程表示如下：

(1)

式中：M∈n×n，為結構系統(tǒng)的質量矩陣;C∈n×n,為結構系統(tǒng)的阻尼矩陣;K∈n×n,為結構系統(tǒng)的剛度矩陣；n×T與n×T，依次為結構系統(tǒng)的位移響應信號、速度響應信號與加速度響應信號的時域采樣矩陣；F(t)∈n×T,為系統(tǒng)的外載荷激勵向量的時域采樣矩陣。

針對阻尼較小的n自由度工程結構，結構的位移響應信號X(t)=[x1(t),x2(t),…,xn(t)]T∈n×T,在模態(tài)坐標中可表示如下：

(2)

式中:矩陣Φ=[φ1,φ2,…,φn]∈n×n,為由結構的n階模態(tài)的振型向量φi(i=1,2,…,n)所構成的模態(tài)振型矩陣；Q(t)=[q1(t),q2(t),…,qn(t)]T∈n×T,為由結構的n階模態(tài)響應qi(t)∈1×T(i=1,2…,n)所組成的模態(tài)響應矩陣。且式(2)中的主振型向量φi需要滿足相互正交，即

(3)

其中mi為第i階模態(tài)質量。并且，qi(t)與qj(t)二者相互獨立，其中i與j為不同的模態(tài)階數(shù)。

工模模態(tài)識別方法的基本思想是通過僅有結構輸出的位移響應信號X(t)，識別出結構的各階模態(tài)振型Φ和模態(tài)響應Q(t)；最后用單自由度識別技術識別出各階的模態(tài)固有頻率f以及模態(tài)阻尼比ξ。

1.2 基于LE算法的工作模態(tài)參數(shù)識別

拉普拉斯特征映射算法是流形學習的熱門算法之一，它的基本思想是從局部的角度去構建數(shù)據(jù)之間的聯(lián)系，這使得降維后的數(shù)據(jù)集依然能保持原來數(shù)據(jù)集的內部結構。LE算法的具體描述如下。

(4)

通過推導式(4)，可得到變換后的LE目標優(yōu)化函數(shù)如下：

mintrace(STLS),

s.t.

STDS=I。

(5)

LS=-DSΛdiag。

(6)

式中：Λdiag為一個對角矩陣，由于D，L均是實對稱矩陣，故它們的轉置與自身相等。因此，對于單獨的si(t)向量，式(6)可寫成以下形式：

Lsi(t)=λDsi(t)。

(7)

LE通過對式(7)做特征值分解，并將其d(d?n)個最小非零特征值[λ1,λ2,…,λi,…,λd]所對應的特征向量{υ1，υ2，…，υi，…，υd}作為降維后的結果輸出，即可達到降維的目的。n個d維向量S(t)={s1(t),s2(t),…,sd(t)}T∈d×T與降維后的結果{υ1，υ2，…，υi，…，υd}一致。

因此，依據(jù)最小二乘廣義逆，存在變換矩陣B=(b1,b2,…,bd)∈n×d，使得式(8)成立，X(t)可分解為式(9)。

B=X(t)ST(t)(S(t)ST(t))-1,

(8)

(9)

1.3 基于LE的三維結構工作模態(tài)參數(shù)識別的物理意義

當使用拉普拉斯特征映射獲取T個S維向量構成的S(t)={s1(t),s2(t),…,sd(t)}T∈d×T(向量si(t)與sj(t)之間滿足相互正交(i≠j))時，由式(2)、式(3)和式(9)，發(fā)現(xiàn)各階模態(tài)響應矩陣Q(t)對應LE得到的低維嵌入S(t)={s1(t),s2(t),…,sd(t)}T∈d×T。模態(tài)振型Φ對應變換矩陣B=(b1,b2,…,bd)∈n×d。因此，模態(tài)參數(shù)識別可通過運用LE分解求得。具體的求解流程如圖1所示。圖1中的流形學習方法可以是線性的主成分分析(PCA)[8]、也可以是非線性的等距離映射(Isomap)[16]。

1.4 基于LE的三維結構工作模態(tài)參數(shù)識別方法

1.4.1 運用最小二乘廣義逆反迭代求解法

針對復雜的三維結構，本文將結構的時域響應信號表示如下：

(10)

在有限元的實驗設計中，復雜三維連續(xù)體結構被離散化為H個檢測點，時間的采樣點總數(shù)為T個。因此，被離散化后的三維連續(xù)結構在時域上的響應在理論上可以近似用式(11)表示。

(11)

式中：S為有限元計算過程中的模態(tài)截斷數(shù)，A(t)S×T為三維系統(tǒng)結構的模態(tài)響應矩陣。LE算法運用最小二乘法廣義逆反代法做三維結構模態(tài)參數(shù)識別的流程圖如圖2所示。

(1)在復雜三維結構的X，Y，Z軸3個方向上的時域數(shù)據(jù)集(XX)H×T(t)，(XY)H×T(t)，(XZ)H×T(t)中選取響應最大的那個方向上的振動數(shù)據(jù)作為識別的基礎數(shù)據(jù)集。將選取的數(shù)據(jù)集記為E,則數(shù)據(jù)集E的維數(shù)為H,樣本點數(shù)為T。

(2)利用LE算法計算數(shù)據(jù)集E降至d維空間上的低維嵌入數(shù)據(jù)A(t)d×T={a1(t),a2(t),…,ai(t),…,ad(t)}T，其中，t=1,2,3,…,T。

(3)用A(t)d×T替代式(11)中的A(t)S×T，用(ΦX)H×d替代式(11)中的(ΦX)H×S，用(ΦY)H×d替代式(11)中的(ΦY)H×S，用(ΦZ)H×d替代式(11)中的(ΦZ)H×S，通過該替換，式(11)可用式(12)表示。

(12)

(4)將步驟(2)中求得的A(t)d×T帶入式(12)中，可依次求出(ΦX)H×d，(ΦY)H×d，(ΦZ)H×d，如式(13)所示:

(13)

1.4.2 直接組裝三維位移響應求解法

最小二乘法廣義逆反代求解法在求解過程中，需要多次求解廣義逆，使得求解精度出現(xiàn)下降，因此本文采用直接組裝三維位移響應求解法。

根據(jù)式(2)，該方法將X，Y，Z三個方向的模態(tài)響應組裝成整個三維結構的整體位移模態(tài)響應。組裝后的公式如式(14)：

(14)

式中:(HThree)3H×T(t)為整個三維結構的整體模態(tài)響應；(XX)H×T(t)，(XY)H×T(t)，(XZ)H×T(t)分別對應X，Y，Z方向的位移響應信號；矩陣Φ3H×S是復雜三維結構的各階模態(tài)主振型向量φi所構成的模態(tài)振型矩陣；QS×T(t)則是結構的各階模態(tài)響應qi(t)所構成的模態(tài)響應矩陣，S是有限元計算過程中的模態(tài)截斷數(shù)。

LE算法運用直接組裝三維位移響應法做三維結構模態(tài)參數(shù)識別的流程圖如圖3所示。

1.5 基于LE的三維結構工作模態(tài)參數(shù)識別與其他方法的比較

本文主要將基于LE的三維結構工作模態(tài)參數(shù)識別方法與基于主成分分析(PCA)[8]、等距離映射算法(Isomap)[16]的三維結構工作模態(tài)參數(shù)識別方法進行比較。正常情況下，算法的時間開銷是不可忽視的問題。除此之外，復雜的三維結構具有一定粘性比例阻尼，因此采集到的數(shù)據(jù)集帶有一些非線性特征。PCA算法是一種線性降維算法，當某階模態(tài)非線性特征過多時，識別復雜三維結構的模態(tài)參數(shù)精度會變得不太理想。LE是一種非線性降維算法，在降維過程中，LE算法能更好地保留數(shù)據(jù)集的特征，從而處理數(shù)據(jù)的非線性特征。因此，LE算法在對復雜三維結構做模態(tài)分析時，某階模態(tài)的識別精度在一定程度上能更優(yōu)于PCA算法。Isomap算法也是非線性降維方法，與LE算法保留局部性質相比，Isomap算法是一種保留全局性質的算法，故這是Isomap時間開銷高的主要原因。在識別精度方面，由于Isomap算法對數(shù)據(jù)的局部結構更為敏感，這會使得其在識別復雜三維結構模態(tài)參數(shù)時，對數(shù)據(jù)的要求更為嚴格，從而造成Isomap的識別精度在一定程度上會低于LE算法。3種方法的優(yōu)缺點說明如表1所示。

表1 基于LE與PCA和Isomap的三維結構工作模態(tài)參數(shù)識別比較

2 仿真驗證結果與分析

2.1 實驗對象和仿真數(shù)據(jù)介紹

為了驗證LE算法能有效識別三維連續(xù)體結構系統(tǒng)的模態(tài)參數(shù)，進行實驗仿真驗證。仿真驗證中主要進行具有復雜三維結構的圓柱殼系統(tǒng)[20-21]研究，模擬圖如圖4所示。

圓柱殼系統(tǒng)有兩端簡支，在其表面放置一定數(shù)量的振動信號傳感器，然后利用設備記錄傳感器X、Y、Z三個方向上的振動響應信號。圓柱殼的參數(shù)設計如下：圓柱殼使用的材料厚度為0.005 m，半徑長度是0.182 5 m，總長0.37 m，同時圓柱殼材料的彈性模量為205 Gpa[8]，泊松比的數(shù)值是0.3，密度是7 850 kg/m3。在進行實驗仿真的過程中，人為地設置阻尼比η為0.03，進行測量。

圓柱殼上的傳感器觀測點布置規(guī)則如下：首先沿著圓柱殼的長度軸向均勻分成38個圈，然后在每個圈中均勻地布置115個觀測點，因此圓柱殼表面共布置H=38×115=4 370個傳感器觀測點。在布置好裝置，通過使用激勵裝置對圓柱殼施加高斯白噪聲激勵的基礎上，再采樣得到三維圓柱殼的振動響應信號。其中，將機器的采樣頻率設置成5 120 Hz，采樣的時間長度設置成1 s，即T=5 120。

振動信號記錄設備采用一套數(shù)采前端LMS SCADS-X Ⅲ，Test lab 9B系統(tǒng)，DELL M65。然后使用LMS Virtual.lab中的有限元方法對傳感器輸出的振動響應信號進行讀取。從而測量到三維圓柱殼結構X、Y、Z三個方向的響應數(shù)據(jù)集合，第358個觀測點的響應信號集合如圖5所示。

2.2 評價指標

為了評價LE算法對三維圓柱體結構系統(tǒng)做模態(tài)參數(shù)識別是否具有有效性，本文將用LE算法識別到復雜三維結構的模態(tài)參數(shù)與利用有限元方法分析識別到的真實模態(tài)參數(shù)進行比較，分析LE算法識別的模態(tài)參數(shù)的有效性。

本文將使用模態(tài)置信度準則(Cmac)來評估LE算法識別到的模態(tài)振型準確度，

(15)

式中:φi為LE算法識別得到復雜三維結構的第i階模態(tài)振型，ψi為利用有限元分析識別復雜三維結構得到的第i階真實模態(tài)振型。由式(15)可知，模態(tài)置信度準則的取值范圍0≤Cmac,i≤1，當Cmac,i的值越接近1時，所識別的第i階的模態(tài)振型精度越高。

2.3 基于LE的三維圓柱體模態(tài)參數(shù)識別結果

本文采用1.4.2節(jié)的直接組裝法識別三維圓柱體的模態(tài)參數(shù)，因為該方法相對最小二乘法廣義逆反代法，對測量數(shù)據(jù)的噪聲不敏感。當圓柱殼阻尼比η=0.03時，為了更直觀地比較，本文將得到的三維模態(tài)振型圖形進行旋轉，得到各階模態(tài)振型的二維形狀圖。本文將識別的模態(tài)振型與有限元分析的真實模態(tài)振型進行比較。有限元分析的真實模態(tài)振型如圖6所示。

圖7是基于LE算法識別的模態(tài)振型三維旋轉成二維圖。LE算法對三維圓柱殼參數(shù)識別出來的固有頻率圖如圖8所示。然后用LE算法識別出的固有頻率與有限元分析計算得到的真實固有頻率進行比較，比較結果如表2所示。最后，根據(jù)模態(tài)置信度準則，求得LE算法的Cmac值，結果如表3所示。

表2 有限元方法與LE算法固有頻率相比較

表3 LE算法置信度準則值

2.4 基于LE的三維圓柱體工作模態(tài)參數(shù)識別結果分析

在1.5節(jié)中，本文采用PCA算法，Isomap算法與LE算法進行了比較，分析了他們各自的優(yōu)缺點。表4和表5給出了利用PCA算法、Isomap算法與LE算法對三維結構模態(tài)參數(shù)識別的比較結果。

表4 PCA算法、Isomap算法與LE算法模態(tài)置信度系數(shù)比較

續(xù)表4

表5 PCA算法、Isomap算法與LE算法固有頻率比較

通過表4可知，與PCA算法相比，LE算法第1、2階的模態(tài)振型識別精度比PCA高，第3、4階差別不大，第7階的識別精度低于PCA。與Isomap算法相比，LE算法除了第1階低于Isomap算法，第2、3、4和7階的識別效果均比Isomap算法好。通過表5可知，LE算法識別第2、3、4階的固有頻率誤差比PCA算法低，第7階固有頻率誤差高于PCA算法。與Isomap算法相比，LE算法識別第3階固有頻率誤差更低，但第7階固有頻率誤差更高。

表6給出了PCA算法、Isomap算法與LE算法的時間開銷，從表6可知，Isomap算法的時間開銷遠高于LE算法，LE算法的時間開銷高于PCA算法。

表6 PCA算法、Isomap算法與LE算法時間開銷比較

3 結束語

本文將拉普拉斯特征映射(LE)算法應用于三維連續(xù)圓柱殼結構系統(tǒng)的模態(tài)參數(shù)識別。通過在圓柱殼結構系統(tǒng)數(shù)據(jù)集上的仿真結果說明，LE算法能有效地識別三維結構系統(tǒng)的工作模態(tài)參數(shù)(模態(tài)振型、模態(tài)響應與固有頻率)。LE算法的識別結果與PCA算法和Isomap算法相比，LE算法能更為均勻地識別出各階模態(tài)振型；且與Isomap相比，LE算法識別速度更快，平均識別精度更高。

但是，LE算法在識別過程中，容易出現(xiàn)某階的模態(tài)缺失或重疊現(xiàn)象。如何解決模態(tài)的缺失或重疊問題，是LE算法模態(tài)參數(shù)識別需要繼續(xù)進一步解決的問題。除此之外，未來的工作中還可以用實際結構進一步驗證該算法的有效性。