任同群,鄭培珍,徐征,王曉東

(1.大連理工大學(xué)精密與特種加工教育部重點實驗室,遼寧 大連116024;2.大連理工大學(xué) 微納米技術(shù)及系統(tǒng)遼寧省重點實驗室,遼寧大連116024)

0 引言

在實際工程中,多數(shù)情況下結(jié)構(gòu)的激勵未知或不易測得,這導(dǎo)致無法準確識別結(jié)構(gòu)的模態(tài)參數(shù),進而限制了基于結(jié)構(gòu)參數(shù)變化的損傷檢測方法的應(yīng)用。相比之下,直接基于振動響應(yīng)信號分析的狀態(tài)檢測方法不影響結(jié)構(gòu)正常運行且費用較低,已成為目前結(jié)構(gòu)健康檢測的研究熱點,也得到了廣泛的工程應(yīng)用,其基本過程包括損傷特征提取與分類[1-4]。

在眾多所提取的損傷特征中,振動傳遞率反映了結(jié)構(gòu)兩測點之間的全部力學(xué)特性,可表示為結(jié)構(gòu)剛度、阻尼和質(zhì)量等物理參數(shù)的函數(shù)。Fauziyah等人[5]利用加速度計陣列測量多泵間的振動傳遞率,實現(xiàn)了多泵的故障診斷。Dong等人[6]基于振動傳遞率對隨機激勵下的螺栓松動監(jiān)測進行了研究,克服了常用的振動方法在發(fā)現(xiàn)螺栓局部松動方面的不足。楊斌[7]采用相空間重構(gòu)方法進行振動傳遞率函數(shù)分析,重構(gòu)振動傳遞率函數(shù)矩陣,并以其奇異值熵作為損傷指標,提高了檢測方法的抗干擾性。Luo等人[8]以加權(quán)傳遞率作為損傷指標評估結(jié)構(gòu)狀態(tài),提高了TF的信噪比。Chen等人[9]利用神經(jīng)網(wǎng)絡(luò),提取基于振動傳遞率的損傷診斷指標,并驗證了其正確性。

在上述研究中,振動傳遞率僅基于單通道信號計算獲得。在工程應(yīng)用中,測試信號是實際振動在傳感器敏感方向的投影。因此,基于單通道信號的方法在處理空間振動問題時容易陷入困難。為解決此問題,通常采用多加速度計或三軸加速度計獲取空間振動信號,進行單通道信號處理及特征提取后再進行特征融合。但是,這種處理方法忽略了信號通道之間的關(guān)聯(lián)性,融合后會帶來失真甚至是虛假的結(jié)果。因此,有必要將多通道信號視為矢量信號,進行同時處理,以保證通道間的關(guān)聯(lián)性。四元數(shù)三個虛部的運算是等價的,因此非常適合表示三維矢量信號。例如,Rehman和Mandic[10]研究了三維信號的經(jīng)驗?zāi)B(tài)分解,并證明三通道信息聯(lián)合處理時會獲得更好的精度。Tong[11]采用四元數(shù)K-LT對故障異步電機進行特征提取,但處理的數(shù)據(jù)僅限于時域四元數(shù)振動序列。Yi和Lv等人[12-13]描述了四通道信號的相關(guān)性,提出一種基于凸優(yōu)化的四元數(shù)奇異譜分析方法。Ren等人[14-15]采用四元數(shù),構(gòu)建三軸信號的振動傳遞率,并驗證了K-LT法和DQWT法,取得了較好的結(jié)果。此外,Ell和Sangwine[16-19]針對四元數(shù)傅立葉變換和相關(guān)分析問題進行了合作研究,相關(guān)研究成果為四元數(shù)頻域分析奠定了理論基礎(chǔ)。

本文提出了一種基于QTJT的狀態(tài)檢測方法,先后采用K-LT和DQWT技術(shù)提取狀態(tài)特征指標向量,然后定義狀態(tài)特征指標向量之間的歐氏距離并將其作為狀態(tài)指標,針對軌道結(jié)構(gòu)扣件松脫及內(nèi)部溫度力的變化,進行了狀態(tài)檢測與識別。通過在無砟軌道試驗平臺上開展實際實驗,驗證該方法的有效性。

1 四元數(shù)聯(lián)合傳遞率相關(guān)理論

1.1 四元數(shù)傳遞率計算

對于單通道信號,振動傳遞率定義為結(jié)構(gòu)兩測試點的時域序列xi(輸出)和xj(輸入)的頻譜比值。它描述了振動(振幅和相位)如何在結(jié)構(gòu)的兩個測試點之間進行傳遞。假設(shè)一個n自由度線性結(jié)構(gòu)受外力f(t)=[f1(t),f2(t),…,fn(t)]的激勵,如果f(t)具有形如白噪聲的均勻譜密度,則fk(t),k=1,2,…,n的傅立葉變換則為一個常數(shù)D。此時,由節(jié)點j到i的振動傳遞率計算公式為

式中:Ai(w),Aj(w)分別為節(jié)點i,j處振動加速度信號的頻譜;H(ω)為n×n的頻響函數(shù)矩陣,H(ω)=-(ω2M+iωC+K)-1;M為n×n的質(zhì)量矩陣;C為n×n的阻尼矩陣;K為n×n的剛度矩陣。如果激振力只作用于一點k,式(1)變?yōu)?/p>

上述兩情形中,激振力僅提供了振動能量,并未參與傳遞率的計算。此時,振動傳遞率僅與激勵位置有關(guān)。根據(jù)式(1)和式(2),振動傳遞率與頻響函數(shù)有關(guān),反映了結(jié)構(gòu)的內(nèi)在特性。因此,振動傳遞率可以用于結(jié)構(gòu)狀態(tài)檢測。

采用三軸加速度計采集獲得結(jié)構(gòu)三個正交方向的振動信號ax(t),ay(t),az(t),可以構(gòu)建純虛四元數(shù)的時域序列aq(t)=[aq1,aq2,…,aqM],其中,aqr=iaxr+jayr+kazr,r∈[1,2,…,M],M為序列維度。此時,傳統(tǒng)的基于單通道信號頻譜之比的振動傳遞率不再適用。此時四元數(shù)傳遞率可定義為

式中:AQ(ω)為四元數(shù)譜;·為四元數(shù)點乘符號;?為共軛運算符號。

Ell提出了四元數(shù)傅立葉變換用于分析二維線性定常系統(tǒng)[16-17],其原始定義可表示為

式(4)定義了二維函數(shù)h(t,τ)從空間域(t,τ)向空間頻域(jω,kv)的轉(zhuǎn)換。四元數(shù)傅立葉變換最早用于彩色圖像的分析,具體形式[18]如下

式中:μ為隨機單位純虛四元數(shù);f(m,n)為M×N彩色圖像矩陣(三顏色通道構(gòu)建純虛四元數(shù));(m,n)和(u,v)分別為空間域和空間頻域的坐標。對于純虛四元數(shù)振動時域序列aq(t),可視其為1×N的彩色圖像,因此,aq(t)的四元數(shù)譜計算公式為

由式(6)可知,AQ(ω)也是四元數(shù)序列,故由頻譜之比定義的四元數(shù)傳遞率也是在頻域內(nèi)的四元數(shù)序列,和普通傳遞率一樣,也包含了全部的結(jié)構(gòu)內(nèi)在特性。

1.2 四元數(shù)傳遞率特性

參考實數(shù)信號的傅立葉變換,四元數(shù)時域序列亦可以由不同頻率分量的四元數(shù)信號的線性組合來表示。類似地,在四元數(shù)傅立葉變換中,可將變換后的模值作為縱坐標得到四元數(shù)幅頻譜。

不失一般性,構(gòu)造如下三通道信號[x(t),y(t),z(t)],三個通道分別由不同頻率的正(余)弦信號組成

對x(t),y(t),z(t)三個單通道信號分別進行獨立的傅立葉變換,得到幅值譜如圖1(a)所示,然后將三個通道信號聯(lián)合成[x(t),y(t),z(t)],構(gòu)建四元數(shù)序列并進行四元數(shù)傅立葉變換,可得到整體的頻域信息,其幅值譜如圖1(b)所示。對比可以發(fā)現(xiàn),聯(lián)合信號的四元數(shù)頻譜包含了全部信號的頻率成分,是三通道信息的無損融合[20]。

圖1 四元數(shù)傅立葉變換Fig.1 Quaternion Fourier transform of three-channel joint signal

當(dāng)振動方向偏離傳感器的敏感方向且僅使用單軸傳感器時,振動幅值會減小,但頻率不變。如果傳感器的安裝方向不能嚴格一致,傳統(tǒng)的基于單通道信號的振動傳遞率將偏離理論值。在實際工程應(yīng)用中,很難保證檢測時傳感器的安裝方向嚴格一致。不僅如此,若振動方向是時變的,將直接導(dǎo)致傳遞率的時變,進而造成結(jié)構(gòu)狀態(tài)的誤診斷。

圖2 (a)描述了任意向量v∈R3繞3D單位向量u旋轉(zhuǎn)2θ的情形。用上節(jié)中的定義將向量u和v擴展為純虛四元數(shù),則旋轉(zhuǎn)之后得到的向量v′計算公式為

式中:euθ為單位四元數(shù)的指數(shù)形式,euθ=cosθ+usinθ。式(8)等價于三軸加速度計的測量坐標系原點不變,坐標系的方向繞向量u旋轉(zhuǎn)2θ,如圖2(b)所示。

圖2 四元數(shù)描述空間向量旋轉(zhuǎn)Fig.2 Quaternion describes vector rotation in space

向量v′的模值計算公式為

上述分析表明當(dāng)采用三軸加速度計時,只要安裝位置固定,無論其安裝姿態(tài)如何變化,均可如實采集獲得空間振動信號。進而可推論,兩測點的三軸加速度計位置固定且不考慮噪聲干擾時,不論安裝姿態(tài)如何,式(3)的四元數(shù)傳遞率幅值均為其理論值。為證明上述推論,基于無砟軌道綜合試驗平臺,獲得兩測點振動信號的四元數(shù)時域序列,通過理論計算將測點A,B的信號分別繞u1,u2旋轉(zhuǎn)π/6和π/3,獲得三組四元數(shù)振動傳遞率,即初原始、僅旋轉(zhuǎn)A信號以及旋轉(zhuǎn)A/B信號的四元數(shù)振動傳遞率,如圖3所示。由圖3可知,當(dāng)傳感器姿態(tài)發(fā)生變化時,四元數(shù)傳遞率的實部與虛部幅值均發(fā)生變化,但整體幅值不變。其中,u1,u2分別為

圖3 信號旋轉(zhuǎn)前后的四元數(shù)傳遞率幅值Fig.3 Amplitude of quaternion transmissibility before and after signal rotation

1.3 對偶樹四元數(shù)小波變換

對偶樹小波變換由離散小波變換(discrete wavelet transform,DWT)和復(fù)小波變換(complex wavelet transform,CWT)拓展而來,基于二維解析信號和希爾伯特變換(Hilbert transform,HT)構(gòu)造而成。與實小波變換相比,對偶樹四元數(shù)小波變換可以獲得多角度、多尺度下的一個幅值和三個相位信息,同時具有良好的時移不變性[21]。

為了推廣到二維,Bülow[22]根據(jù)四元數(shù)傅里葉變換給出了四元數(shù)解析信號的定義,給定一個實數(shù)二維信號f(x,y),其四元數(shù)解析信號定義為

式中:fHi1(x,y)和fHi2(x,y)分別為沿x和y軸的希爾伯特變換;fHi(x,y)為完整的希爾伯特變換。

二維離散小波變換由一維離散小波變換沿兩個正交基的張量積構(gòu)成,可表示為尺度函數(shù)φ(x)φ(y)、水平方向子帶函數(shù)φ(x)ψ(y)、豎直方向子帶函數(shù)ψ(x)φ(y)和對角線方向子帶函數(shù)ψ(x)ψ(y)的組合[23]。相似地,四元數(shù)小波同樣具有一個尺度函數(shù)和三個小波子帶函數(shù),即將二維離散小波變換與三個實小波(由一維希爾伯特變換得到)組合構(gòu)成一個四元數(shù),從而得到二維解析小波。以對角子帶ψ(x)ψ(y)為例,其二維希爾伯特變換計算公式為

式中:(ψh,ψg)為特定濾波器h和g對應(yīng)的小波函數(shù),它遵循{ψh,ψg=HT(ψh)}的原則。式(12)中的每個分量均可由一維對偶樹復(fù)小波的組合計算。因此,對角線子帶四元數(shù)小波計算公式為

同理可得水平、垂直子帶四元數(shù)小波以及四元數(shù)尺度函數(shù)[24]分別為

對偶樹小波變換的計算細節(jié)可參考文獻[25]。

2 基于四元數(shù)聯(lián)合傳遞率的狀態(tài)檢測方法

2.1 卡洛變換(K-LT)法

K-LT是一種基于對象統(tǒng)計特征的最優(yōu)正交變換。它可以消除統(tǒng)計相關(guān)性,但保留分類信息。將狀態(tài)已知條件下測量得到的四元數(shù)傳遞率模值向量作為訓(xùn)練樣本,并構(gòu)建為狀態(tài)矩陣SQ=[s1,s2,…,sN]∈RM×N,其中M為測量四元數(shù)傳遞率的長度,即參與運算的有效譜線數(shù)目,N為測量的傳遞率數(shù)目。狀態(tài)矩陣SQ的協(xié)方差矩陣CQ∈RM×N可以分解為

式中:W=[w1,w2,…,wN]∈RM×N為特征向量矩陣張成的特征子空間;Σ=diag(λ1,λ2,…,λN)∈RN×N為特征值矩陣,且λ1≥λ2≥…≥λN。對任意訓(xùn)練四元數(shù)傳遞率模值向量sj,j=1,2,…,N,將其投影到特征子空間可獲得新的向量sj′=WT·sj。此時,s′j可作為狀態(tài)索引向量。將未知狀態(tài)下測得的四元數(shù)傳遞率模值向量以同樣方式投影到特征子空間,獲得投影向量s′t,然后計算s′j和s′t的歐氏距離,稱取得最小歐式距離的投影向量為匹配向量,即同一狀態(tài)。

2.2 對偶樹四元數(shù)小波變換(DQWT)法

如圖4所示,將同一狀態(tài)下4個測試四元數(shù)傳遞率模值向量疊加,構(gòu)建狀態(tài)特征圖。圖4表明同一狀態(tài)下的不同樣本有很大的相似性而又有著局部的差異性。

圖4 四元數(shù)傳遞率灰度圖Fig.4 Grey-scale map of quaternion transmissibility

對狀態(tài)特征圖進行對偶樹小波變換,其系數(shù)可構(gòu)成一個矩陣F,即

F的第一列包含圖像的低頻信息,第2~4列分別代表圖像在水平方向、垂直方向和對角線方向的高頻信息。

對偶樹小波變換的系數(shù)包含了一個幅值信息和三個相位信息。F的每一列可以構(gòu)成式(17)的四元數(shù),計算其每個四元數(shù)的幅值和相位,并重新組合得到幅值-相位矩陣

獲得每種已知狀態(tài)下的多組四元數(shù)傳遞率,構(gòu)建狀態(tài)特征圖作為訓(xùn)練樣本,對其進行對偶樹四元數(shù)小波變換得到低頻幅值相位向量AFi,同樣地,構(gòu)建未知狀態(tài)下的狀態(tài)特征圖,計算其低頻幅值相位向量AFj,計算AFi和AFj之間的最小歐式距離并將其作為狀態(tài)識別指標,與上節(jié)方法類似。取低頻幅值相位信息,既綜合考慮了二維傳遞率圖像的幅值相位信息,又減小了單傳遞率的測量不確定性對檢測結(jié)果的影響,且降低了原始數(shù)據(jù)的大小。

3 狀態(tài)檢測實驗

實驗對象為總長20 m的CRTS-I型“無砟軌道綜合試驗平臺”,鋼軌兩端鎖緊,中間由配套的WJ-7型扣件鎖緊。實驗時在軌底隔跨兩組扣件安裝兩個三軸加速度傳感器PCB 356A33,采用力錘YD-5T在實驗區(qū)域左側(cè)某處的軌頂施加激勵,敲擊方向垂直于軌頂并保持激勵位置固定,以10 kHz的采樣頻率采集振動信號,采集設(shè)備為DHDAS-5902動態(tài)信號采集分析系統(tǒng)。實驗現(xiàn)場及傳感器布置如圖5所示。

圖5 實驗現(xiàn)場及傳感器布置Fig.5 Experiment site and sensor location

首先,利用鋼軌扣件的松脫與否模擬不同狀態(tài),根據(jù)扣件松脫順序設(shè)置五種不同工況,分別為:扣件全部鎖緊,扣件1松脫,扣件1和2松脫,扣件1,2和3松脫以及4個扣件全部松脫。每種工況下,在左側(cè)激勵點用力錘激勵40次。根據(jù)前述四元數(shù)傳遞率函數(shù)的計算方法,計算得到每種工況下的40組傳遞率函數(shù)。圖6為工況1下的全部40次四元數(shù)傳遞率,取200~1700 Hz部分傳遞率函數(shù),可見其重合性很好[26]。一方面,說明了傳遞率計算的正確性;另一方面,也驗證了四元數(shù)傳遞率函數(shù)僅由激勵位置和結(jié)構(gòu)自身特性確定,而與激勵大小無關(guān)的結(jié)論。

圖6 工況1下40組原始四元數(shù)傳遞率函數(shù)幅值譜[26]Fig.6 Amplitude spectrums of 40 sets of original quaternion transmissibility under situation 1[26]

分別求取五種工況下的40組四元數(shù)傳遞率函數(shù)平均值,如圖7所示。在選取頻段內(nèi),傳遞率譜線大體重合,僅在部分頻率區(qū)間存在較明顯差別。這是因為相較實驗區(qū)域,試驗平臺在較長一段鋼軌上的扣件是鎖緊的,僅松脫局部相鄰扣件對整體結(jié)構(gòu)的影響較小。

圖7 五種工況下平均四元數(shù)傳遞率函數(shù)幅值譜Fig.7 Average amplitude spectrum of quaternion transmissibility under five situations

采用2.1中基于K-LT的方法進行狀態(tài)識別,分別取五種工況的前30組傳遞率函數(shù)作為訓(xùn)練樣本,后10組作為測試樣本,用來驗證該方法的有效性。檢測結(jié)果如表1所示,其中,測試QTJT(i)下的第一列數(shù)字表示與之取得最小歐氏距離的訓(xùn)練樣本列序號,正確識別時第j測試組下的QTJT(i)對應(yīng)序號應(yīng)為30×(j-1)+1~30j;第二列表示狀態(tài)識別結(jié)果,‘√’表示識別正確,‘×’表示識別錯誤。由表1可知,該方法取得了98%的狀態(tài)識別正確率。

表1 基于K-LT方法扣件松脫狀態(tài)檢測結(jié)果[26]Tab.1 Detection results based on K-LT method

基于上述數(shù)據(jù),采用2.2中DQWT方法進行狀態(tài)識別處理,分解層次取2,每種狀態(tài)下計算得到的40組QTJTs模值矩陣V∈R152×40,每4組QTJT作為一個樣本,每種狀態(tài)得到10組樣本,其中前6組作為訓(xùn)練樣本,后4組作為測試樣本。識別結(jié)果如表2所示,狀態(tài)i(i=1,2,3,4,5)正確的匹配范圍應(yīng)為6(i-1)+1~6i。由表2可知,該方法狀態(tài)識別正確率達到了100%。

表2 基于DQWT方法扣件松脫狀態(tài)檢測結(jié)果[26]Tab.2 Detection results based on DQWT method

采用LG-900型液壓鋼軌拉伸機拉伸鋼軌,模擬鋼軌內(nèi)部縱向溫度力的變化,分別設(shè)置拉應(yīng)力:0,10,20,30,40,50 MPa。首先,松脫全部扣件以釋放鋼軌內(nèi)部縱向溫度力,穩(wěn)定后拉伸機分別以上述設(shè)置應(yīng)力動作進行鋼軌拉伸,而后再將全部扣件鎖緊,縮緊力矩為150 N·m。同樣地,每次拉伸之后按照前面的實驗方法,依次松脫相鄰兩組扣件,得到30組狀態(tài),每種狀態(tài)采集20組振動信號,計算得到20組四元數(shù)傳遞率。在每種狀態(tài)下,以其前12組QTJT分別構(gòu)建3個狀態(tài)特征圖作為訓(xùn)練樣本,以后8組QTJT分別構(gòu)建2個狀態(tài)特征圖作為測試樣本,表3為該實驗的狀態(tài)組織安排情況?；贒QWT方法進行狀態(tài)識別處理,每種狀態(tài)的匹配范圍應(yīng)為3(i-1)+1~3i,所有樣本的識別結(jié)果如表4所示,該方法狀態(tài)識別正確率達到了100%。

表3 溫度力、扣件松脫綜合實驗結(jié)構(gòu)狀態(tài)安排Tab.3 Structure state arrangement for comprehensive experiment of temperature force and fastener looseness

表4 的實驗結(jié)論表明,該方法可實現(xiàn)分辨力優(yōu)于10 MPa的鋼軌溫度力變化識別。基于鐵路軌道相關(guān)標準,鋼軌的相關(guān)參數(shù)如下:截面積為7.745×10-3m2,線膨脹系數(shù)α=0.0118 mm/(m·℃),彈性模量E=2.1×105 MPa。則軌溫有1℃變化時所產(chǎn)生的內(nèi)部縱向力計算公式為

表4 基于DQWT溫度力、扣件松脫綜合鋼軌狀態(tài)檢測結(jié)果[26]Tab.4 Detection results based on DQWT method for the comprehensive experiment

LG-900型液壓鋼軌拉伸機的油缸截面積為S=3.352×10-3m2,則其10 MPa的拉應(yīng)力對應(yīng)10 MPa×S=33520 N的拉力。進而,10 MPa的拉應(yīng)力相對應(yīng)的鋼軌溫度變化為33520 N/19192.11 N×1℃=1.75℃。因此,本方法可實現(xiàn)分辨力優(yōu)于10 MPa的鋼軌溫度力變化識別,相當(dāng)于鎖定狀態(tài)下1.75℃的軌溫變化。

4 結(jié)論

將傳統(tǒng)的基于單通道信號的振動傳遞率推廣到多通道信號的四元數(shù)振動傳遞率,并結(jié)合K-LT技術(shù)和DQWT,給出了基于四元數(shù)聯(lián)合傳遞率的的結(jié)構(gòu)狀態(tài)檢測方法。通過在無砟軌道試驗平臺上的實際試驗,驗證了該方法的有效性,狀態(tài)識別結(jié)果正確率高于95%。上述方法提供了一種識別狀態(tài)變化是否存在的解決方案,未來可以通過增加測試點和依次識別相鄰四元數(shù)傳遞率,進一步發(fā)展為同時識別狀態(tài)變化的存在和位置。更重要的是,該方法保留了振動傳遞率的理論值,對傳感器的安裝方向沒有嚴格的要求,在實際應(yīng)用中具有明顯優(yōu)勢。特別是基于DQWT的方法,綜合考慮了多個被測傳遞率,降低了單次測量不確定性的影響。本文的研究為推動直接基于振動響應(yīng)信號分析的狀態(tài)檢測方法發(fā)展起到了積極作用,具有重要技術(shù)應(yīng)用價值。