李廣敏,陳江偉,紀(jì)曉輪,張旭野,王星星
(江蘇國電南自海吉科技有限公司,江蘇南京 210009)
在工業(yè)領(lǐng)域中大量使用旋轉(zhuǎn)機(jī)械,如風(fēng)機(jī)、水泵、發(fā)電機(jī)等,這些設(shè)備一旦出現(xiàn)故障,往往會(huì)造成經(jīng)濟(jì)損失甚至安全事故。為防止旋轉(zhuǎn)設(shè)備嚴(yán)重故障的發(fā)生,工業(yè)領(lǐng)域中越來越多的采用狀態(tài)監(jiān)測及故障診斷技術(shù)[1]。
目前工業(yè)診斷中通常采用設(shè)置閾值或數(shù)值的變化趨勢實(shí)現(xiàn)故障預(yù)警,再通過大量的故障診斷算法實(shí)現(xiàn)故障識(shí)別。但是不同種類和大小的旋轉(zhuǎn)設(shè)備需要研究不同的振動(dòng)診斷算法,這使得故障診斷的通用性很差,降低了設(shè)備故障識(shí)別的準(zhǔn)確率,往往當(dāng)確定設(shè)備存在故障后,需要通過人力來實(shí)現(xiàn)診斷[2]。
通過對旋轉(zhuǎn)設(shè)備常見故障的研究發(fā)現(xiàn),相同類型故障的信號(hào)頻譜之間存在相似性,通過將設(shè)備信號(hào)的頻譜與典型故障頻譜進(jìn)行比對,容易識(shí)別出設(shè)備的故障類型,并且此工作可由計(jì)算機(jī)完成。本文主要從數(shù)學(xué)的角度,借助于統(tǒng)計(jì)手段,研究任意曲線幾何形態(tài)相似性的定義、度量與應(yīng)用問題。
定義:設(shè)T為相似閾值,a為兩種頻譜的相似度,且a和T在0到1閉區(qū)間取值,當(dāng)a>T時(shí),認(rèn)為兩種故障頻譜曲線相似,T為經(jīng)驗(yàn)值。
設(shè)圖1曲線為f(x),圖2曲線為g(x),如果圖1與圖2曲線幾何形態(tài)相似,則必滿足f(x)=k×g(x)(k>0)。
圖1 f(x)曲線圖Fig.1 Curve of f(x)
圖2 g(x)曲線圖Fig.2 Curve of g(x)
可以得出曲線f(x)與g(x)相似的必要條件為:f(x)=k×g(x)(k>0)。
下面論述f(x)=k×g(x)(k>0)是曲線f(x)與g(x)相似的充分條件。引用歐幾里得幾何里關(guān)于三角形相似的說法:存在兩個(gè)三角形,其中一個(gè)三角形通過放大或縮小n倍后能夠與另外一個(gè)三角形完全重合,那么這兩個(gè)三角形相似。由該定理擴(kuò)展到曲線相似可得:存在兩條曲線f(x)與g(x),當(dāng)其中一條曲線經(jīng)過放大或縮小后能夠完全與另外一條曲線重合,那么這兩條曲線相似。
f(x)=k×g(x)(k>0)說明g(x)經(jīng)過k倍變化后能與f(x)完全重合,所以f(x)=k×g(x)(k>0)是曲線f(x)與曲線g(x)相似的充分條件。并且定義當(dāng)曲線f(x)與曲線g(x)相似時(shí),相似度為1。
當(dāng)f(x)=k×g(x)(k>0)在x的定義域范圍內(nèi)不能完全成立時(shí),那么曲線f(x)與曲線g(x)不完全相似,引入相似度數(shù)學(xué)量表示兩曲線的相似程度[3]。
由f(x)=k×g(x)得k=f(x)/g(x)
對于定義域范圍內(nèi)的任意x成立,設(shè)集合X={x0,x1,x2,x3,…,xn}為 x 在定義域范圍內(nèi)的取樣序列,則可得出:
當(dāng)兩曲線不完全相似時(shí)此連等式將不成立,由此推導(dǎo)相似度計(jì)算方法,圖1與圖2為兩曲線,X軸為頻率,Y軸為幅值。將兩曲線等間隔的分割,將曲線看成由多個(gè)矩形組成,設(shè)R1(n)表示圖1中自左向右第n個(gè)矩形高度,R2(n)表示圖2中自左向右第n個(gè)矩形高度,兩曲線分割的矩形數(shù)量相等為 m,設(shè)[4]:
幾何中三角形相似定理指出,兩個(gè)三角形對應(yīng)三條邊成比例,即可得出兩個(gè)三角形相似,由此以曲線被分割的小矩形類比三角形的邊給出曲線相似度計(jì)算的公式,設(shè):
相似度a=1-[R1(1)/R2(1)+R1(2)/R2(2)+R1(3)/R2(3)+… +R1(m)/R2(m)]/m。通過比較a與T的大小判斷兩波形是否近似相似[5]。
圖3為函數(shù)f(t)=sin(2×pi×15×t)+2×sin(2×pi×40×t)的時(shí)域波形和經(jīng)過FFT變換后頻譜波形,圖4為函數(shù)g(t)=0.8×sin(2×pi×40×t)+0.3×sin(2×pi×15×t)的時(shí)域波形和經(jīng)過FFT變換后的頻譜波形,圖中采樣點(diǎn)m=128,振幅單位為數(shù)學(xué)約化單位1(下同)。通過函數(shù)表達(dá)式和頻譜圖可以定性看出兩者是相似的。
圖3 f(t)時(shí)域與頻域波形Fig.3 f(t)in frequency domain and time domain waveform
設(shè)F(f)為f(t)的頻譜函數(shù),G(f)為g(t)的頻譜函數(shù),f為頻率,f在定義域[0,200]內(nèi)等間隔取樣128點(diǎn),計(jì)算可得集合A={F(f0),F(xiàn)(f1),F(xiàn)(f2),…,F(xiàn)(f128)},B={G(f0),G(f1),G(f2),…,G(f128)}
設(shè) D={A0/B0,A1/B1,A2/B2,…,A128/B128},則相似度 a=1-sum(D)/128=0.9069。可以說明f(t)與g(t)的頻譜有90%的相似度,該計(jì)算方法可以很好的給出兩曲線相似度的度量方法。
圖4 f(x)曲線圖Fig.4 g(t)in frequency domain and time domain waveform
圖5 k(t)時(shí)域與頻域波形Fig.5 k(t)in frequency domain and time domain waveform
對f(t)作出修改,使k(t)=2×sin(2×pi×15×t)+1×sin(2×pi×40×t),圖(5)為該函數(shù)的時(shí)域和頻域波形。按以上方法計(jì)算修改后的f(t)與g(t)的頻譜相似度為0.4200,說明兩者相似度很低。通過圖(5)與圖(4)的對比可以看出兩函數(shù)在時(shí)域與頻域上的相似度也是很低的,說明此相似度計(jì)算方法可以很好地反應(yīng)曲線間相似程度。
由計(jì)算公式可以看出,兩曲線的相似度值受計(jì)算時(shí)m的取值影響,表1列出m為128、256、512、1024、1536時(shí)函數(shù) f(t)=sin(2×pi×15×t)+2×sin(2×pi×40×t)與函數(shù) g(t)=0.8×sin(2×pi×40×t)+0.3×sin(2×pi×15×t)頻域的相似度值。
通過對表1數(shù)據(jù)的分析,m的取值變化對相似度值影響不大,建議工程運(yùn)用中采用多種m值求相似度平均值的方法。
表1 相似度與采樣點(diǎn)數(shù)的關(guān)系Table 1 The relationship between similarity and the number of sampling points
本文基于幾何相似和統(tǒng)計(jì)學(xué)原理,給出了任意時(shí)域波形的頻譜相似度的定義與度量方法,通過應(yīng)用算例計(jì)算所得數(shù)據(jù),驗(yàn)證了該算法在時(shí)域波形頻譜相似度計(jì)算的準(zhǔn)確性。工程應(yīng)用中結(jié)合旋轉(zhuǎn)設(shè)備故障頻譜庫,可以有效的識(shí)別設(shè)備故障類型。
[1]沈慶根.設(shè)備故障診斷[M].北京:化學(xué)工業(yè)出版社,2006.
[2]楊國安.滾動(dòng)軸承故障診斷實(shí)用技術(shù)[M].北京:中國石化出版社,2012.
[3]文海玉.應(yīng)用數(shù)理統(tǒng)計(jì)[M].哈爾濱:哈爾濱工業(yè)大學(xué)出版社,2012.
[4]胡慶平.BCI-代數(shù)[M].西安:陜西科技出版社,1987.
[5]Pawlak Z.Rough set[J].International Journal of Computer and System Sciences,1993,46:44-54.