肖 峻，高東俊，奚修豪

(武漢理工大學(xué) 機(jī)電工程學(xué)院，湖北 武漢 430070)

鹽鹵是制作豆腐的關(guān)鍵成分，也是工業(yè)中提取生產(chǎn)鹵素的重要原料。鹵水具有強(qiáng)腐蝕性，且在運(yùn)輸過程中溫度較高，使得鹽鹵泄漏不僅會對環(huán)境造成嚴(yán)重污染，還會給人們帶來巨大的經(jīng)濟(jì)損失[1]，故需要采用有效的方法對鹵水管道泄漏進(jìn)行及時檢測和準(zhǔn)確定位，從而使經(jīng)濟(jì)損失和環(huán)境污染最小化。隨著新興技術(shù)的涌入，如壓力變送器技術(shù)[2]、光纖技術(shù)[3-4]、聲波法[5]、小波變換技術(shù)[6]、經(jīng)驗(yàn)?zāi)B(tài)分解[7-8]等，負(fù)壓波法的定位精度與泄漏檢測準(zhǔn)確率得到提升[9]。為此，筆者采用改進(jìn)的小波閾值消噪算法對采集的壓力信號消噪，基于極值的方法獲取壓力突變時間標(biāo)簽點(diǎn)，并通過優(yōu)化的定位算法進(jìn)行定位計算，以期提高鹵水管道泄漏檢測與定位效果。

1 負(fù)壓波檢測原理

由于負(fù)壓波頻率較低，在管道內(nèi)的反射衰減速度較慢，所以負(fù)壓波可沿管道兩端傳播很遠(yuǎn)的距離[10]。當(dāng)負(fù)壓波到達(dá)管道兩端時，會出現(xiàn)一段時間壓力陡降，這時可通過負(fù)壓波到達(dá)管道兩端的時間差和負(fù)壓波在管道中的傳播速度來確定泄漏點(diǎn)的位置，故可在管道兩端分別布置壓力傳感器，采集負(fù)壓波信號。負(fù)壓波檢測原理如圖1所示。

負(fù)壓波泄露點(diǎn)定位包括負(fù)壓波到達(dá)管道兩端傳感器的時間差和負(fù)壓波波速兩個關(guān)鍵因素，泄漏點(diǎn)定位公式為：

(1)

式中：x為泄漏點(diǎn)到管道首端的距離；L為管道長度；α為負(fù)壓波傳播速度；V為管道內(nèi)流體流速；Δt為泄漏所產(chǎn)生的負(fù)壓波傳到管道首、末兩端的時間差。

2 負(fù)壓波波速計算

2.1 負(fù)壓波波速通用計算公式

負(fù)壓波的傳播速度隨著管道物理參數(shù)、傳輸介質(zhì)物理性質(zhì)及外界物理參數(shù)的變化而具有不同的表現(xiàn)形式，此處考慮純液中負(fù)壓波的傳播速度，根據(jù)質(zhì)量守恒、液體的連續(xù)性和液體壓縮與管壁膨脹原理，負(fù)壓波傳播速度可表示為[11]：

(2)

式中：α為負(fù)壓波傳播速度；ρ為管道液體密度；K為管道液體體積模量；E為管道材料的彈性模量；D為管道的平均內(nèi)徑；δ為管道的管壁厚度；ψ為管道修正系數(shù)。

2.2 負(fù)壓波波速公式修正

在實(shí)際工程中，溫度是影響負(fù)壓波傳播速度的首要物理因素，液體密度和液體體積模量都會隨溫度變化而變化。由于鹵水的粘度較小，在管道內(nèi)的流動性好，所以溫度在鹵水管道中傳播比較均勻。通過采集現(xiàn)場的鹵水，將鹵水加熱到不同溫度并對其密度進(jìn)行測量，得到不同溫度下的鹵水密度，并對其進(jìn)行線性擬合。其中，溫度在65～74 ℃范圍內(nèi)鹵水密度與溫度的關(guān)系如圖2所示。

圖2 鹵水密度與溫度的關(guān)系

進(jìn)一步得到鹵水密度與距離的表達(dá)式為：

(3)

式中：x為泄漏點(diǎn)到管道首端的距離；τ1、τ2分別為管道首端與末端的溫度測量值；L為管道長度。

現(xiàn)場測得鹵水礦化度為10 g/L，根據(jù)文獻(xiàn)[12]中同礦化度鹵水在不同溫度下體積模量的測量數(shù)據(jù)，得到體積模量與距離的表達(dá)式為：

(4)

從而得到修正后的負(fù)壓波傳播速度：

(5)

3 壓力信號消噪

3.1 消噪算法優(yōu)化

小波軟閾值消噪算法處理后的信號具有很好的連續(xù)性，但處理后的小波系數(shù)與原小波系數(shù)之間存在一定偏差，會給重構(gòu)后的信號帶來誤差[13-14]。為此，采用指數(shù)型閾值函數(shù)對閾值函數(shù)進(jìn)行優(yōu)化，其表達(dá)式如式(6)所示。

(6)

以Bumps加噪信號為例進(jìn)行消噪仿真分析，軟閾值消噪算法與指數(shù)型閾值函數(shù)消噪算法相關(guān)性能指標(biāo)對比如表1所示，消噪效果對比圖如圖3所示，結(jié)合表1和圖3可知，指數(shù)型閾值函數(shù)消噪算法整體降噪效果較好，且在較關(guān)鍵的最大突變段起始點(diǎn)與終止點(diǎn)處波形特征體現(xiàn)更好，便于后期讀取時間標(biāo)簽點(diǎn)。

表1 軟閾值消噪算法與指數(shù)型閾值函數(shù)消噪算法性能對比

圖3 消噪對比圖

3.2 時間標(biāo)簽點(diǎn)的讀取

在已消噪的壓力信號基礎(chǔ)上讀取負(fù)壓波中幅值變化最大的陡降段的首、末兩個極值點(diǎn)，并將其作為時間標(biāo)簽點(diǎn)，原因在于負(fù)壓波在傳播至管道首、末兩端壓力傳感器過程中會衰減，但衰減主要體現(xiàn)在幅值上，其頻率變化較小，整個波形的時間標(biāo)簽點(diǎn)相對固定。為了得到更準(zhǔn)確的泄漏時間差，可將管道首、末兩端壓力傳感器所檢測到的壓力最大突變段起始點(diǎn)與結(jié)束點(diǎn)作為時間標(biāo)簽點(diǎn)，根據(jù)實(shí)驗(yàn)數(shù)據(jù)進(jìn)行采樣點(diǎn)的截取，確定時間標(biāo)簽點(diǎn)，具體如圖4所示。

圖4 時間標(biāo)簽點(diǎn)讀取示意圖

圖4中，a、b分別表示管道首、末兩端圧力曲線的最大突變起始點(diǎn)，c、d分別表示管道首、末兩端圧力曲線的最大突變結(jié)束點(diǎn)。則泄漏時間差的計算公式為：

(7)

式中：n1、n2、n3、n4分別為a、b、c、d4個點(diǎn)所對應(yīng)的采樣點(diǎn)序列數(shù)；f為系統(tǒng)采樣率。

4 泄漏點(diǎn)的定位計算

4.1 泄漏點(diǎn)定位算法

根據(jù)式(5)，某次泄漏所產(chǎn)生的負(fù)壓波傳播到管道首、末兩端所用時間分別為：

(8)

式中：t1和t2分別為負(fù)壓波傳播到首、末兩端的時間；xlack為泄漏點(diǎn)坐標(biāo)。

筆者采用Romberg算法計算數(shù)值積分，以t1的數(shù)值積分為例，首先將積分區(qū)間[0,xlack]分成n等份。設(shè)x0=0,xn=xlack,則有：

xk=kh

(9)

(10)

將步長減半，得到積分域2n等份后的近似值為：

(11)

于是可以得到遞推公式為：

(12)

故由截斷誤差表達(dá)式可得：

(13)

式中：RTn(α)、RT2n(α)分別為選取步長h、h/2計算得出的截斷誤差；I(α)為負(fù)壓波由泄漏點(diǎn)傳遞至管道首端的實(shí)際時間；Tn、T2n分別為選取步長h、h/2計算得出的負(fù)壓波由泄漏點(diǎn)傳遞至管道首端的時間。

雖然在管路中負(fù)壓波傳播速度隨位置不同而發(fā)生變化，但波速關(guān)于位置的二階導(dǎo)數(shù)相對變化不大，故α″(ξ1)≈α″(ξ2)。由式(13)可得：

(14)

(15)

由于α(x)關(guān)于x的四階導(dǎo)數(shù)和六階導(dǎo)數(shù)變化不大，可通過Simpson值、Cotes值到Romberg值的迭代計算，加速收斂求解過程，3個值的計算方法如式(16)所示。若精度要求為ε，則當(dāng)|R2n-Rn|≤ε時，可達(dá)到精度要求。

(16)

4.2 算法優(yōu)化

雖然通過Romberg算法對上、下變積分的計算更準(zhǔn)確，但其搜索范圍較廣泛，故首先壓縮搜索范圍，然后進(jìn)行區(qū)間二分搜索。

4.2.1 壓縮搜索范圍

假定管路中波速為常數(shù)值，將首、末兩端計算所得波速α(0)、α(L)代入式(1)進(jìn)行計算，得到泄漏點(diǎn)分別為xf、xe。由于α(x)為單調(diào)遞減函數(shù)，所以實(shí)際泄漏點(diǎn)應(yīng)在(xf,xe)區(qū)間內(nèi)。

4.2.2 區(qū)間二分搜索

在(xf,xe)區(qū)間內(nèi)，根據(jù)1/2比例對搜索范圍進(jìn)行不斷壓縮。每一次代入的位置計算點(diǎn)為上一次搜索所確定區(qū)間的中點(diǎn)，并由新的搜索點(diǎn)來判斷下一次的搜索區(qū)間。區(qū)間二分搜索示意圖如圖5所示。將xi(xi∈(xf,xe),i=1，2，…，m)代入到Romberg算法中計算近似時間值t1和t2，求得時間差為Δti，進(jìn)一步求得Δti與實(shí)測時間差Δt的差值。

圖5 區(qū)間二分搜索示意圖

以點(diǎn)x2為例，取x1=(xf+xe)/2，可得到y(tǒng)(x1)的值，然后比較y(x1)、y(xf)的大小。若y(xf)×y(xe)<0且y(xe)×y(x1)>0，則舍去(x1,xe)段，再對(xf,x1)區(qū)間進(jìn)行二分搜索，如圖5實(shí)線所示，取x2=(xf+x1)/2為下一次二分搜索點(diǎn)。若y(xf)×y(x1)>0且y(xe)×y(x1)>0，則舍去(xf,x1)段，再對(x1,xe)區(qū)間進(jìn)行二分搜索，如圖5虛線所示，下一次區(qū)間二分搜索點(diǎn)為x2=(xe+x1)/2。按照上述方法不斷進(jìn)行二分搜索，直到搜索區(qū)間的寬度小于或等于設(shè)定閾值δ，取最小搜索區(qū)間的中點(diǎn)為最終定位點(diǎn)。

4.3 定位算法實(shí)現(xiàn)與性能分析

運(yùn)用Matlab對優(yōu)化后的Romberg算法編程，泄漏點(diǎn)定位計算流程圖如圖6所示。假設(shè)管道長度L=20 km、時間差Δt=15 s，波速表達(dá)式相同，定位精度分別為1%和5%的情況，進(jìn)行100次仿真實(shí)驗(yàn)，得到定位計算耗時情況如表2所示，可看出兩種定位精度的算法耗時基本相同，對泄漏點(diǎn)的定位計算實(shí)時性影響不大，故選取定位精度為1%進(jìn)行定位計算。

表2 定位耗時比較

圖6 泄漏點(diǎn)定位計算流程圖

5 實(shí)驗(yàn)與結(jié)果

5.1 實(shí)驗(yàn)檢測裝置

為了驗(yàn)證上述方法的可行性，在室外環(huán)境下搭建模擬鹵水管道泄漏檢測與定位試驗(yàn)平臺，選取管徑D為20 mm、壁厚δ為1 mm的薄壁管道，實(shí)驗(yàn)中所使用液體為工業(yè)鹽溶解液，管道壓力范圍為0～0.8 MPa，通過開關(guān)球閥來模擬泄漏。管道結(jié)構(gòu)設(shè)置為U型，設(shè)定兩壓力傳感器間管道總長度為33.25 m，兩球閥的安裝位置與首端壓力傳感器的距離分別為12.00 m和21.25 m。實(shí)驗(yàn)采用STM32單片機(jī)和Matlab軟件平臺對壓力信號進(jìn)行采集和分析，具體如圖7所示。

圖7 室外模擬實(shí)驗(yàn)平臺

5.2 實(shí)驗(yàn)結(jié)果與分析

實(shí)驗(yàn)選取泄漏點(diǎn)為12.00 m和21.25 m的數(shù)據(jù)進(jìn)行分析和計算，現(xiàn)場采集的壓力信號與消噪情況如圖8所示，泄漏定位結(jié)果如表3所示。結(jié)合圖8和表3可以看出，壓力信號消噪效果良好，泄漏點(diǎn)的定位精度較高。

圖8 壓力信號消噪情況

表3 泄漏定位結(jié)果

6 結(jié)論

利用負(fù)壓波法可以對長直輸鹵管道進(jìn)行泄漏點(diǎn)檢測和定位。筆者通過對鹵水密度與溫度關(guān)系的標(biāo)定，修正了負(fù)壓波波速計算公式，減少了溫度對波速的影響。采用小波閾值消噪算法對采集數(shù)據(jù)進(jìn)行濾波，有效降低了噪聲干擾，有助于拾取壓力突變點(diǎn)。選取Romberg積分算法進(jìn)行定位計算，大幅提高了計算準(zhǔn)確度。實(shí)驗(yàn)結(jié)果表明，上述研究方法能有效提高鹵水管道泄漏檢測和定位效果。