韓文翔, 章蘭珠
(華東理工大學 機械與動力工程學院,上海 200237)
管道作為運輸介質的載體,具有經濟、安全、不間斷性、運輸量大等多重優勢。在提供便利的同時,由于腐蝕、老化和人為破壞等因素的影響,一旦發生泄漏事故,后果頗為嚴重,損壞國民經濟甚至威脅生命安全[1]。針對管道泄漏的檢測已經有示蹤法、負壓波法、聲學檢測法等較成熟的檢測方法。其中聲發射(acoustic emission,AE)技術作為一種具有在線動態監測特點的新型檢測技術,可以及時發現管道的早期缺陷,提高檢測效率,具有良好的工程實際應用價值[2]。
目前在管道泄漏領域的研究較多集中在單點泄漏場景,利用時差定位法、神經網絡等方法可有效解決管道泄漏的定位問題。而真實工況下往往會發生多點泄漏的現象,多點泄漏由于泄漏信號之間的相互干擾給研究帶來了一定的難度,目前還沒有形成成熟的管道多點泄漏定位方法。高華[3-4]等提出一種改進差分進化算法的管道泄漏定位方法,將采集到的信號經過小波消噪后并根據信號幅值建立改進差分進化算法的目標優化函數,通過對其求解獲得泄漏源位置,但是該方法未考慮彎管、變徑等工況對信號幅值的干擾。而聲發射信號可以避免復雜工況對信號特性的影響,管道在發生泄漏時由于內外壓差導致內部介質與外部環境之間的摩擦從而形成高頻應力波[5],這種機械波組成的聲發射信號與語音信號同為非平穩信號,可借鑒語音信號處理中的盲辨識法應用于管道多點泄漏的檢測。盲辨識問題實質上是盲源分離中混合矩陣的估計問題,根據信號源數量與傳感器數量之間的關系可分為超定、適定和欠定3種類型[6]。目前該領域的研究熱點主要在欠定盲源分離,而超定和適定兩種情況下的盲源分離技術相對成熟,本文的研究重點在于利用信號處理技術提取有效的泄漏特征信號并結合盲源分離算法實現管道的多點泄漏源定位,因此選用基于快速獨立成分分析(fast independent component analysis,Fast ICA)的適定盲源分離算法。
在盲源分離的應用中,往往將噪聲分量當做獨立成分分離,考慮到去噪效果和盲源分離后排序未知等問題的限制,本文先對采集信號進行去噪預處理。通過對變分模態分解(variational mode decomposition,VMD)分解后有效模態選取依據的研究提出一種結合互譜的信號去噪算法,并通過試驗驗證了該去噪方法的有效性和多點泄漏應用的可靠性。
聲發射技術是近年來得到快速發展的一種動態無損檢測新技術,該技術可以在不破壞被檢結構的情況下對構件進行全天在線監控,并能對結構中出現的疲勞及失效做出早期預警,具有實時監測、動態監測、覆蓋面廣等優點[7]。管道發生泄漏時,釋放的氣體與管道裂口通過摩擦所產生的能量以應力波的形式釋放,此過程稱為聲發射現象,通過傳感器對聲發射信號的采集和分析可實現管道泄漏的檢測。但是由于聲發射信號的本身特性及環境因素的影響,需要對采集信號進行適當的分析和處理才能進一步精確確定泄漏點的位置。
在泄漏源定位方面,常采用時差定位法[8],其原理如下
(1)
式中:x為被檢測管道泄漏點的位置,即泄漏點到上游傳感器距離,m;Δt為泄漏信號到達兩個傳感器的時間差,s;v為泄漏信號在管道中的傳播速度,m/s。
由于兩傳感器之間的距離是已知的,因此只需知道Δt便可得出泄漏源的位置。互相關函數是計算兩個信號相關性的常用方法,根據互相關時域分析圖中最大相關峰值進而計算得到上,下游有效泄漏信號到達兩個傳感器的時差Δt。
管道運行過程中環境噪聲的干擾會對定位結果造成顯著影響,因此有必要先對采集到的信號進行去噪處理,目前信號去噪領域已具很多成熟的算法,例如小波分解、小波包分解等信號處理算法。雖然小波分解在工業中的應用已經有了長遠的發展和進步,但是小波基函數的選擇和分解層數的確定始終是擺在分析人員面前的問題[9]。為此Huang等[10]提出一種完全基于數據驅動的自適應分解方法——經驗模態分解(empirical mode decomposition,EMD),它將信號分解成有限個從高頻到低頻的固有模態函數(intrinsic mode function,IMF) 和殘余分量之和, 能夠很好地分析處理非線性、非平穩信號。由于信號中存在間斷信號和脈沖干擾等異常事件的干擾,使得EMD 存在模態混疊的問題,這會降低EMD 的降噪效果?;诖?,研究人員又陸續發明了集合經驗模態分解(ensemble empirical mode decomposition,EEMD),VMD等信號處理算法[11]。
互譜可由互相關函數經傅里葉變換求得從而在頻域內描述兩信號s1(k)和s2(k)之間的相關性[12]。對s1(k)和s2(k)的相關函數作傅里葉變換得互譜公式如下
Ss1s2(f)=αSss(f)ej2πfτ=αSss(f)ejθs1s2(f)
(2)
式中,Sss(f)為采集信號的自譜。由式(2)可得時間延遲為τ,則相位差為
θs1s2(f)=2πfτ
(3)
在采集信號的頻帶內,θs1s2(f)與頻率成線性關系,斜率為2πτ?;プV的計算公式如下
(4)
式中:S1(f)和S2(f)分別為s1(k)和s2(k)的離散傅里葉變換系數;*為復共軛。
互譜估計常采用Welch平均周期圖法表示為
Gs1s2(f)=Gs1s2(f)+jQs1s2(f)
(5)
式中,Gs1s2和Qs1s2(f)分別為共譜密度和正交譜,可得相位譜如下
(6)
將s1(k)和s2(k)的互譜相位譜作差可得互譜相位差譜,差譜值越接近水平表明頻率相關性越高。
利用Welch平均周期圖法可得S1(f)和S2(f)的自譜Gx1x2(f)和Gx2x2(f),則s1(k)和s2(k)的幅值平方相干函數為
(7)
相干函數可從頻域的角度反應兩段信號的相關程度,頻帶的相干函數值越大表明該頻段的信噪比越高,因此可根據互譜相位差譜和相干函數估計泄漏特征頻帶。
試驗采用北京軟島科技研發的聲發射儀器和實驗室搭建的管道裝置,在泄漏源上下游兩端各放置一個聲發射傳感器,選用2.5 M的采樣率,計算兩段信號的互譜相位差譜和相關函數如圖1所示。

圖1 s1(t)和s2(t)的相干函數及互譜相位差譜Fig.1 Coherence function and cross-spectral phase difference spectrum of s1(t) and s2(t)
由圖1(a)可知,相關性較高的頻段主要分布在0~500 kHz;由圖1(b)可0~500 kHz內未發現明顯水平差譜頻段,分析其原因可能是泄漏特征信號與噪聲信號相互重疊從而影響了互譜效果。為此提出一種互譜結合變分模態分解的去噪算法。
VMD通過迭代搜尋變分模型最優解來確定每個分量的頻率中心及帶寬,相比于EMD,VMD具有堅實的理論基礎,有效地彌補了 EMD和EEMD等不足[13]。采集信號經過VMD后得到若干個IMF(FIM)分量和殘余分量,公式如下
(8)
式中:k為VMD分解后的IMF分量個數;γ(t)為分解后的殘余分量。
2.2.1 分解層數的確定
相比于EMD自適應分解信號的特點,VMD的分解層數K由人為輸入,如果K過小,低頻段將被劃分為高頻段,給模式函數的求解帶來了困難。但是如果K值過大,會導致過度分解,計算時間過長。因此,分解層數的確定是VMD在應用過程中必須考慮的問題[14]。
由于采樣率選用2.5 MHz,選用經驗值K=10將采集信號進行VMD分解后肘頻圖,如圖2所示。

圖2 VMD分解結果時頻圖Fig.2 Time-frequency diagram of VMD decomposition result
由圖2(b)可知,采集信號從低頻到高頻依次分解,并且只要分解層數選取適當可有效避免模態混疊等問題。由于高頻部分的分解毫無意義,不僅消耗算力并且容易造成過度分解。因此可以通過互譜大致定位有效頻段并結合濾波實現信號的去噪預處理,通過巴特沃斯濾波器濾出其中的高相似度信號,再對濾波后的信號進行VMD分解。
將信號進行互譜分析并設計低通濾波器濾出0~500 kHz的頻域圖,如圖3所示。

圖3 濾波后頻域圖Fig.3 Frequency domain graph after filtering
可通過是否發生過度分解為界限確定分解層數,對濾波后的信號進行4層分解和5層分解,得出的頻譜圖如圖4所示。

圖4 VMD分解4層和5層頻域圖Fig.4 Result diagram of VMD decomposition layer 4 and layer 5
由圖4可知,當分解至5層時第4個IMF和第5個IMF出現了過度分解現象,因此將分解層數定位4層。
2.2.2 有效模態選取
基于時差定位的管道泄漏定位試驗中,主要的泄漏特征信號往往包含在VMD分解后的一個或多個模態中,因此針對有效模態的選取進行了深入研究。
目前已有的研究中,往往通過相關系數、能量值和信息熵[15]等熵值判斷有效模態的選取,而以上方法都有一定的局限性。由于高頻噪聲的干擾,通過相關系數和能量值選取有效模態往往導致結果誤差較大。而后者依據的理論是聲發射作為一種動態檢測方法,泄漏特征成分越多,對應信號的隨機性越強,從而信息熵、樣本熵等熵值較大,該方法具有一定的適用性,但是在泄漏率低的情況下仍然無法辨識出有效模態。
將采集信號經過VMD分解成10個模態后,分別計算每個模態的能量值,相關系數和信息熵結果如圖5所示。
由圖5可知,能量值最高的模態為IMF1,信息熵最大的模態是IMF6,相關系數最大的模態是IMF4??梢娀谛畔㈧睾拖嚓P系數選擇的的IMF頻段并不是泄漏特征信號所在的主要頻段。為保證該試驗的科學性,選取0.1 L/min(小泄漏),0.5 L/min(中泄漏),1 L/min(大泄漏)3種泄漏率下各10組試驗,按經驗值將信號分解至10層后,分別計算各個模態的能量值,相關系數和信息熵,記錄每個模態得到最大值的次數,得到的結果如圖6~圖8所示

圖6 0.1 L/min 泄漏率下各指標對應模態選取系數Fig.6 Modal selection coefficient of each index under 0.1 L/min leakage rate

圖7 0.5 L/min 泄漏率下各指標對應模態選取系數Fig.7 Modal selection coefficient of each index under 0.5 L/min leakage rate

圖8 1.0 L/min 泄漏率下各指標對應模態選取次數Fig.8 Number of modal selection corresponding to each index under 1.0 L/min leakage rate
由圖6~圖8可知,基于單一指標的情況下容易出現高頻模態被選用的情況從而導致定位結果存在較大的誤差。為此提出一種結合相關系數、信息熵和能量值的綜合指標,其中相關系數的選用為避免噪聲的干擾采用基于互譜濾波后的信號與分解后的IMF進行互相關計算。將3個指標分別歸一化后相加,取值最大的IMF為包含最多泄漏特征信號的IMF。
設能量值為a,信息熵為b,相關系數為c,將每個指標對應的數值進行歸一化后再將不同指標的數值相加即得到綜合指標xzonghe如下
(9)
信號經過VMD并通過綜合指標選取IMF的結果分布如圖9所示。

圖9 不同泄漏率下綜合指標選取次數Fig.9 Selection times of comprehensive indicators under different leakage rates
由圖9可知,采用綜合指標后選取的模態基本上集中在低頻模態,結果與理論泄漏特征信號所在頻段相符合,所以采用綜合指標選取有效的IMF比單一指標具有更高的可靠性。
通過以上信號處理方式保證了去噪后的信號包含了盡可能多的泄漏特征信號,進而在此基礎上通過FastICA算法將不同泄漏源產生的泄漏信號分離出來。
FastICA其目標就是計算矩陣A的逆矩陣B,使其作用在觀測信號x(t)上,從而可以分離解耦出具有非高斯性的接近原始信號s(t)的源信號[16],基本原理如圖10所示。

圖10 ICA模型Fig.10 ICA model
FastICA算法有基于峭度、基于似然最大、基于負熵最大等形式,這里采用基于負熵最大的FastICA算法,該算法以負熵最大作為一個搜尋方向,可以實現順序地提取獨立源,可以有效地將不動點迭代所帶來的優良算法特性與負熵所帶來的更好統計特性結合起來,使得收斂更加快速、穩健。
綜合以上分析,本文提出一種互譜結合VMD的信號處理方法,并選用綜合指標最大值作為選取有效IMF的依據,在管道同一側對去噪后的信號進行盲源分離計算得到兩個包含不同泄漏源特征的信號,再根據相關系數匹配出管道兩側由同一泄漏源產生的信號,最后根據時差定位得到定位結果,完整的泄漏定位方法如圖11所示。

圖11 多點泄漏定位算法Fig.11 Multi-point leak location algorithm
為證明所提出方法的有效性,研究以兩泄漏點為例,利用北京軟島科技公司研發的聲發射采集設備對實驗室搭建的管道泄漏模型展開研究,在管道上、下游分別放置兩個聲發射傳感器,通過轉動泄漏閥模擬泄漏源,設置兩端傳感器距離為1 300 mm,以左端傳感器為基準,兩泄漏源位置分別為300 mm和1 000 mm,泄漏率大小分別為1 L/min和0.5 L/min,設置管道內壓為0. 3 MPa。試驗方案如圖12所示。

圖12 多點泄漏定位試驗方案(mm)Fig.12 Experimental scheme of multi-point leak location(mm)
試驗裝置如圖13所示。

圖13 管道泄漏模型Fig.13 Pipeline leakage model
以上游的兩個傳感器為例,兩傳感器采集的混合信號s1和s2頻域圖如下圖14、圖15所示。

圖14 混合信號s1頻域圖Fig.14 Frequency domain diagram of mixed signal s1

圖15 混合信號s2頻域圖Fig.15 Frequency domain diagram of mixed signal s2
將信號進行互譜結合VMD去噪并根據綜合指標選取的IMFs1和IMFs2頻域圖如圖16、圖17所示。

圖16 去噪信號IMFs1頻域圖Fig.16 IMFs1 frequency domain diagram of denoising signal

圖17 去噪信號IMFs2頻域圖Fig.17 IMFs2 frequency domain diagram of denoising signal
對IMF1和IMF2進行FastICA計算后得到兩個獨立的信號源q1和q2如圖18、圖19所示。

圖18 信號源q1頻域圖Fig.18 Frequency domain diagram of signal source q1

圖19 信號源q2頻域圖Fig.19 Frequency domain diagram of signal source q2
同理計算管道下游采集的兩路信號s3和s4,經過去噪處理和FastICA計算后得到信號源q3和q4,對q1,q2和q3,q4兩兩計算相關系數值如表1所示。

表1 信號源相關系數計算
由表1可知,(q1,q3)和(q2,q4)的相關系數較高,分別計算其相關延遲時間并得到最終的定位結果。在該工況下進行10組試驗得到的試驗結果如圖20、圖21所示。

圖20 泄漏源位置為300 mm的定位結果Fig.20 Locate the leakage source at 300 mm

圖21 泄漏源位置為700 mm的定位結果Fig.21 Locate the leakage source at 700 mm
為衡量該方法的可靠性,對試驗結果進行數據分析,箱形圖可以用來觀察數據整體的分布情況,利用中位數,25%分位數,75%分位數,上邊界,下邊界等統計量來描述數據的整體分布情況?;谠囼灲Y果繪制箱型如圖22、圖23所示。

圖22 泄漏源位置為300 mm的定位統計結果Fig.22 Locate the leakage source at 300 mm

圖23 泄漏源位置為700 mm的定位統計結果Fig.23 Location statistics for leak source 700 mm
由圖22、圖23可知,泄漏源為300 mm的定位結果較為密集得分布在300 mm左右,誤差率在2%~15%,泄漏源為700 mm的定位結果存在一個偏離分布的異常值,可將其舍棄,總體誤差率在1%~10%?;谝陨辖Y果分析可得,該多點泄漏定位方法具有較高的可靠性和有效性。
本文通過對信號進行變分模態分解后模態選取方法的研究,提出一種結合多指標的模態選取方法,并根據互譜和變分模態分解各自的局限性將兩種算法結合應用,利用該算法對采集信號進行去噪預處理,將去噪后的信號通過FastICA實現源信號的分離,進而通過時差定位實現管道的多點泄漏定位。
最后,通過實驗室的管道泄漏模型開展了多點泄漏定位試驗,通過試驗結果的分析驗證了該方法的有效性,為工業實況提供了一種管道多點泄漏的理論思路。