李赟杰,周 領(lǐng),歐傳奇,劉德有
(1.河海大學(xué)水利水電學(xué)院,江蘇 南京 210098; 2.國(guó)際小水電中心,浙江 杭州 310002)
管道漏損廣泛存在于流域輸水、城市供水等工程中,泄漏不僅造成巨額的經(jīng)濟(jì)損失,而且可能造成次生災(zāi)害,如影響水質(zhì)、損壞建筑基礎(chǔ)等,因此尋求準(zhǔn)確有效的泄漏檢測(cè)方法是一個(gè)重要的研究課題。目前,泄漏檢測(cè)方法主要有目視觀測(cè)法、探地雷達(dá)檢漏法、示蹤劑檢漏法、振動(dòng)聲學(xué)檢漏法、光纖傳感技術(shù)法、瞬態(tài)檢漏法[1]以及穩(wěn)態(tài)水力模型法等,其中,瞬態(tài)檢漏法是泄漏檢測(cè)的主要手段[2-3]。泄漏主要通過(guò)增加阻尼率和在產(chǎn)生的壓力中創(chuàng)建新的反射信號(hào)兩個(gè)方面影響瞬態(tài)信號(hào),識(shí)別和量化這些影響是瞬態(tài)泄漏檢測(cè)技術(shù)的核心[4]。
Pudar等[1]在1992年第一次提出基于壓力測(cè)量的檢漏程序和逆問(wèn)題的解決方法,即已知壓力和流量求管道粗糙度和泄漏情況。誤差分析表明系統(tǒng)摩擦系數(shù)估計(jì)的誤差帶來(lái)的影響比水頭測(cè)量誤差的影響大得多,摩擦系數(shù)的精確性是逆瞬態(tài)法的關(guān)鍵[5-7]。目前,為了提高泄漏檢測(cè)的精確度,除了較大泄漏問(wèn)題,水力瞬變模型中一般均需要考慮動(dòng)態(tài)摩阻的影響[2,8-10]。Zielke[11]推導(dǎo)出適用于管道瞬變層流的卷積加權(quán)模型,該模型考慮了管道瞬時(shí)平均流速和加權(quán)歷史流速變化,被認(rèn)為是精確動(dòng)態(tài)摩阻模型;Vardy等[12-14]在Zielke[11]的基礎(chǔ)上建立了光滑管和粗糙管內(nèi)流動(dòng)的瞬態(tài)摩擦加權(quán)函數(shù)模型,將卷積加權(quán)模型推廣到阻力平方區(qū);Brunone等[15]在Daily等[16]研究基礎(chǔ)上考慮對(duì)流加速度,引入動(dòng)態(tài)摩阻系數(shù),提出了適用于層流和紊流的動(dòng)態(tài)摩阻模型,該模型能很好地平衡運(yùn)算量和計(jì)算精度的問(wèn)題;Bergant等[17]和Vitkovsky等[18]對(duì)Brunone模型分別進(jìn)行了公式改進(jìn)和近似求解,使Brunone模型更加簡(jiǎn)潔、適用。郭新蕾等[3]根據(jù)實(shí)際工程需要研究了Brunone模型的離散網(wǎng)格和求解方法,將其應(yīng)用于管道泄漏模擬,壓力波形對(duì)比表明該模型比傳統(tǒng)的準(zhǔn)恒定泄漏摩阻模型更加精確。然而,加權(quán)函數(shù)類動(dòng)態(tài)摩阻模型很少用于管道泄漏檢測(cè)的研究。
本文將瞬態(tài)加速類(Brunone模型)和加權(quán)函數(shù)類(Zielke模型和Vardy模型)兩類動(dòng)態(tài)摩阻模型用于水力瞬變泄漏建模和模擬,并將計(jì)算結(jié)果與已有的試驗(yàn)數(shù)據(jù)和模擬結(jié)果進(jìn)行對(duì)比,采用3種模型進(jìn)行了算例計(jì)算分析,通過(guò)對(duì)比不同瞬變流模型泄漏過(guò)程中的壓力和流量變化情況,分析泄漏畸變點(diǎn)突變?cè)?、Brunone模型計(jì)算誤差和泄漏中部分參數(shù)的敏感性。
常見的動(dòng)態(tài)摩阻模型主要有以下兩類[10,19]:①基于加權(quán)函數(shù)考慮歷史速度和歷史加速度對(duì)當(dāng)前流態(tài)影響的加權(quán)函數(shù)類模型,如Zielke模型和Vardy模型;②基于瞬時(shí)加速度的模型,如Brunone模型及其改進(jìn)型。本文將Zielke模型、Vardy模型和Brunone模型用于管道泄漏檢測(cè)的水力瞬態(tài)建模。
Zielke[11]在考慮管道瞬時(shí)平均流速和加權(quán)歷史流速變化的情況下,推導(dǎo)出了適用于層流的動(dòng)態(tài)摩阻項(xiàng)JU:
(1)
Vardy等[12]建立了光滑管內(nèi)流動(dòng)的瞬態(tài)摩擦加權(quán)函數(shù)模型,并假定湍流黏度在均勻流速核心周圍的剪切層內(nèi)呈線性變化,對(duì)于高雷諾數(shù)的流動(dòng),給出了水力光滑管道動(dòng)態(tài)摩阻模型:
(2)
(3)
式中:x為沿管道中心線方向的距離;a為波速;k3為Brunone動(dòng)態(tài)摩阻系數(shù),取值的適當(dāng)性對(duì)模型精確程度有很大影響,k3的確定方法有兩種:一是由試驗(yàn)確定;二是采用Vardy等[12]所給的系數(shù)公式:
(4)
Bergant等[17]發(fā)現(xiàn)Brunone模型動(dòng)態(tài)摩阻項(xiàng)計(jì)算公式具有局限性:對(duì)于某些加速或減速條件,該公式對(duì)流項(xiàng)符號(hào)錯(cuò)誤,因此他們?cè)谄浠A(chǔ)上提出了改進(jìn)模型:
(5)
Vitkovsky等[18]隨后給出了該模型的顯性和隱性近似解法。
管道泄漏的數(shù)學(xué)模型由控制方程、動(dòng)態(tài)摩阻和邊界條件3部分組成。本文采用特征線法對(duì)管道泄漏模型進(jìn)行求解。動(dòng)態(tài)摩阻的求解一般包含隱式和顯式求解兩種格式[3]:隱式求解格式中,將穩(wěn)態(tài)摩阻項(xiàng)和動(dòng)態(tài)摩阻項(xiàng)共同作為管道流中摩擦引起的水頭損失代入雙曲型微分方程,這種方法改變了特征線網(wǎng)格,運(yùn)算復(fù)雜且求解時(shí)間長(zhǎng),優(yōu)點(diǎn)是更為精確;顯式求解格式中,將動(dòng)態(tài)摩阻作為附加項(xiàng)單獨(dú)計(jì)算,可以較好地代入各種模型,模型簡(jiǎn)單且計(jì)算時(shí)間短。綜合考慮模型精度和運(yùn)算效率,本文采用第二種方法。
管道瞬變流的連續(xù)性方程和動(dòng)量守恒方程[20-22]如下:
(6)
(7)
式中:H為測(cè)壓管水頭;J為管道水體運(yùn)動(dòng)摩擦引起的單位長(zhǎng)度的水頭損失;f為Darcy-Weisbach摩擦系數(shù),其值大小依賴于雷諾數(shù)和管道相對(duì)粗糙度。
如圖1所示,將線性雙曲型微分方程(6)(7)變換為常微分方程,采用一階有限差分近似,分別沿著正負(fù)特征線進(jìn)行積分,得到如下C+、C-特征線方程:
(8)
(9)
式中:vP、HP分別為在t時(shí)刻節(jié)點(diǎn)P處的流速和水頭;vR、HR分別為在t-Δt時(shí)刻節(jié)點(diǎn)R處的流速和水頭;vS、HS分別為在t-Δt時(shí)刻節(jié)點(diǎn)S處的流速和水頭;JR、JS分別為在t-Δt時(shí)刻節(jié)點(diǎn)R、S處的水頭損失;Δt為計(jì)算時(shí)間間隔。
圖1 特征線積分網(wǎng)格
上游端為已知水位的水庫(kù),下游為一閥門,泄漏孔口邊界條件如圖2所示。
圖2 泄漏孔邊界條件
假定泄漏發(fā)生在第i個(gè)節(jié)點(diǎn)處,泄漏孔口外部為大氣壓,正負(fù)特征線測(cè)壓管水頭H1、H2與泄漏處測(cè)壓管水頭相等。由孔口出流公式得泄漏孔口的流量為
(10)
式中:Cd為泄漏孔口的流量系數(shù),可由試驗(yàn)確定具體值,取值范圍一般在0.60~0.61之間[23];Ag為泄漏孔口面積;ΔHL為泄漏點(diǎn)的壓力水頭,ΔHL=HL-zL,其中HL為管道泄漏處的測(cè)壓管水頭,zL為管道泄漏點(diǎn)高程。
由連續(xù)性方程可得:
Q1=Q2+QL
(11)
式中:Q1、Q2分別為泄漏節(jié)點(diǎn)處的正向特征線和負(fù)向特征線對(duì)應(yīng)的流量。將式(10)(11)與式(8)(9)聯(lián)立求解,可得流量Q1、Q2、QL和測(cè)壓管水頭H1、H2、HL。
不考慮泄漏情況,將本文動(dòng)態(tài)摩阻模型(Zielke模型、Vardy模型和Brunone模型)的模擬結(jié)果與Bergant等[17]的實(shí)測(cè)數(shù)據(jù)進(jìn)行對(duì)比。試驗(yàn)裝置基本參數(shù)如下:管道長(zhǎng)L=37.23 m,直徑D=0.022 1 m,上游水箱恒定水頭h=32 m,波速a=1 319 m/s,管道初始流速v0=0.3 m/s,雷諾數(shù)Re=5 600。閥門瞬間關(guān)閉,模擬的閥門處壓力(用測(cè)壓管水頭Hv表示)如圖3所示。由圖3可知,與不考慮動(dòng)態(tài)摩阻的傳統(tǒng)摩阻模型相比,動(dòng)態(tài)摩阻模型能精準(zhǔn)描述瞬變管道壓力幅值和衰減,隨著周期數(shù)增加,動(dòng)態(tài)摩阻模型的模擬值與實(shí)測(cè)值仍能吻合較好,而傳統(tǒng)摩阻模型模擬的壓力波在兩個(gè)周期后誤差逐漸增大,相位偏移加大。
圖3 不同摩阻模型模擬的閥門處測(cè)壓管水頭與實(shí)測(cè)值對(duì)比
考慮管道泄漏情況,將本文模擬結(jié)果與Lee等[24]模擬的結(jié)果進(jìn)行對(duì)比,其模擬參數(shù)為:L=2 000 m,波速a=1 000 m/s,D=0.5 m,h=50 m,初始管道流量Q0=0.02 m3/s,Re=50 000,f=0.02,泄漏孔口距上游水庫(kù)300 m,初始泄漏量QL=0.2Q0。不考慮動(dòng)態(tài)摩阻的影響,模擬結(jié)果見圖4。從圖4可以看出,兩次模擬的壓力幅值和衰減幾乎完全一致,從而證明了本文模型的正確性。
圖4 泄漏情況下閥門處壓力對(duì)比
如圖5所示,選擇水庫(kù)-管道-閥門(RPV)模型進(jìn)行管道泄漏瞬變模擬。假定管道為一水平直管,模型參數(shù)為:L=2 000 m,D=0.2 m,上游水庫(kù)水位恒定為60 m,a=1 000 m/s,Q0=0.02 m3/s,Re=85 000,f=0.03。
圖5 泄漏管道概化物理模型
假設(shè)泄漏孔外界為大氣壓強(qiáng),泄漏位置(距上游水庫(kù)的距離)XL=1 000 m,QL=0.1Q0。節(jié)點(diǎn)數(shù)N=100,t=0時(shí)末端閥門瞬間關(guān)閉,計(jì)算時(shí)間為32 s,分別模擬了該泄漏情況下的閥門處壓力Hv和管道入口處流量Q,結(jié)果如圖6所示。
圖6 有泄漏時(shí)不同摩阻模型模擬的閥門處壓力和管道入口處流量
由圖6可知,管道出現(xiàn)泄漏時(shí),管中壓力和流量會(huì)出現(xiàn)畸變。相比于流量,壓力畸變更加容易觀察,因此使用壓力畸變信號(hào)能更好地探尋泄漏點(diǎn)。
以泄漏時(shí)閥門處的壓力變化為例,閥門關(guān)閉前,管道中水體處于恒定流狀態(tài),水庫(kù)水位保持恒定,壓強(qiáng)為p0,當(dāng)閥門突然關(guān)閉時(shí),緊靠閥門處的水體流速瞬間變?yōu)榱悖瑝簭?qiáng)立刻升高Δp,上游水體碰到靜止不動(dòng)的水體時(shí)就像碰到完全關(guān)閉的閥門一樣,類似狀態(tài)以波速a向上游傳播,這種波稱為主壓力波,主壓力波對(duì)管道壓力和流量變化起主要作用。當(dāng)主壓力波到達(dá)泄漏點(diǎn)時(shí),泄漏處的壓力發(fā)生突然變化,在泄漏點(diǎn)形成向泄漏兩側(cè)傳播的次壓力波,次壓力波也會(huì)影響管道壓力變化。圖7給出了管道泄漏時(shí)該工況下不同時(shí)刻的壓力波運(yùn)動(dòng)圖形。
由圖7(a)可知t=2 s時(shí)閥門處壓力畸變的原因?yàn)椋褐鲏毫Σ◤拈y門向上游傳播,經(jīng)過(guò)泄漏點(diǎn)時(shí)形成新的次壓力波。主壓力波從閥門傳遞至泄漏點(diǎn)前引起水體壓力增大,當(dāng)這種增壓波遇到泄漏點(diǎn)(假定為大氣壓)時(shí),由于壓差作用會(huì)向閥門處反射一個(gè)引起水體壓力減小的減壓信號(hào),在t=2 s時(shí)到達(dá)閥門,引起閥門處壓力突降。
圖7 泄漏時(shí)不同時(shí)刻管道的壓力波運(yùn)動(dòng)
由圖7(b)可知t=4 s時(shí)閥門處壓力驟降的原因?yàn)椋涸鰤翰◤拈y門向水庫(kù)傳播,引起逐段的水體壓力增大,傳播至泄漏點(diǎn)時(shí),不僅有向閥門處的減壓波,也會(huì)產(chǎn)生向水庫(kù)方向的次壓波,這種波也為減壓波。由于主壓力波起主導(dǎo)作用,主增壓波和次減壓波疊加仍然會(huì)使泄漏點(diǎn)上游的管道壓力增加,傳播至水庫(kù)時(shí),由于壓力差Δp的作用,水庫(kù)向閥門反射一個(gè)使下游壓力降低的波,在t= 3 s時(shí)又到達(dá)泄漏處產(chǎn)生新的次波,此時(shí)仍有從泄漏點(diǎn)反射的其他次波,同樣的疊加作用使得泄漏點(diǎn)下游管道壓力降低,t= 4 s時(shí)到達(dá)閥門導(dǎo)致閥門處壓力驟降。
Zielke模型主要適用于層流,也可用于雷諾數(shù)不高(Re<106)的湍流;Vardy水力光滑管道動(dòng)態(tài)摩阻模型[12]適合于湍流情況,Vardy改進(jìn)模型[13]給出了雷諾數(shù)在2×103~108情況下更為精確的系數(shù)計(jì)算公式;Brunone模型在層流和紊流時(shí)均可適用。Zielke模型和Vardy模型通常被認(rèn)為是精確動(dòng)態(tài)摩阻模型,Brunone模型有簡(jiǎn)化,存在一定的簡(jiǎn)化誤差[25]。在XL=1 000 m,QL=0.1Q0工況下得到壓力如圖8所示,其中圖8(a)為3種動(dòng)態(tài)摩阻模型模擬的閥門處實(shí)際壓力對(duì)比,圖8(b)是以Vardy模型為標(biāo)準(zhǔn),3種模型模擬的閥門處的壓力分別減去Vardy模型的模擬值所得到的壓力對(duì)比曲線。
圖8 不同動(dòng)態(tài)摩阻模型模擬的泄漏時(shí)閥門處壓力對(duì)比曲線
由圖8(a)可知,Brunone模型在前半個(gè)周期內(nèi)所模擬的壓力幅值和衰減與Zielke模型、Vardy模型基本一致,隨著時(shí)間增長(zhǎng),因泄漏導(dǎo)致的壓力畸變值雖趨于一致,但其整體壓力衰減值小于Zielke模型和Vardy模型,兩個(gè)周期后的壓力曲線與Zielke模型和Vardy模型產(chǎn)生明顯差異。圖8(b)可以更好地印證圖8(a)中所觀察到的現(xiàn)象,Brunone模型在t=4 s時(shí)開始與Zielke模型和Vardy模型產(chǎn)生明顯偏離,偏離值隨時(shí)間增長(zhǎng)而增大,并每隔半個(gè)周期出現(xiàn)一次較大的異變點(diǎn),Zielke模型和Vardy模型在圖8(a)(b)中均無(wú)明顯差異??赡艿脑蚴呛?jiǎn)化的Brunone模型有一定的模型誤差,Brunone模型的系數(shù)k3對(duì)模型的精確程度有較大影響,需要現(xiàn)場(chǎng)率定。
圖9 不同泄漏量影響下的管道壓力、流量變化
以Vardy模型為例,考慮泄漏模型在不同泄漏量和不同泄漏位置情況下的管道壓力、流量變化。模型參數(shù)不變,泄漏位置距上游水庫(kù)1 000 m,初始泄漏量QL分別取0.01Q0、0.05Q0、0.1Q0和0.2Q0,圖9比較了4種工況下閥門處壓力Hv、泄漏處壓力HL、管道入口處流量Q和泄漏處流量Q1。由圖9可見,初始泄漏量越大,閥門起始?jí)毫υ叫?,這可由水錘基本公式驗(yàn)證。隨著泄漏量增加,壓力和流量阻尼增大,收斂速度加快,壓力畸變點(diǎn)隨著泄漏量增加變得越來(lái)越清晰,畸變量與初始泄漏量呈正相關(guān),與時(shí)間呈負(fù)相關(guān)。
為了分析不同泄漏位置情況下的管道壓力、流量變化,固定QL=0.1Q0不變,泄漏位置XL分別取距上游水庫(kù)500 m、1 000 m、1 500 m,模擬結(jié)果如圖10所示。由圖10可知,不同泄漏位置會(huì)導(dǎo)致壓力和流量畸變點(diǎn)位置變化,以閥門壓力第一周期為例,泄漏點(diǎn)越靠近閥門,畸變?cè)皆纾儠r(shí)間基本與泄漏距閥門距離成正比,據(jù)此可根據(jù)畸變時(shí)間反推泄漏位置,其他壓力和流量畸變情況類似,可由壓力波運(yùn)動(dòng)情況推得;泄漏量相同,泄漏位置不同時(shí)壓力和流量的幅值不同,阻尼不同,泄漏位置越靠近閥門,壓力和流量曲線幅值越小,衰減越快。
圖10 不同泄漏位置影響下的管道壓力、流量變化
a. 動(dòng)態(tài)摩阻模型比傳統(tǒng)準(zhǔn)恒定摩阻模型更加精確,能更好地模擬管道泄漏時(shí)管中各點(diǎn)的壓力和流量變化。與Zielke和Vardy動(dòng)態(tài)摩阻模型相比,Brunone模型存在一定的簡(jiǎn)化誤差,模型系數(shù)須現(xiàn)場(chǎng)率定。
b. 與流量畸變信號(hào)相比,在管道中用壓力畸變信號(hào)探尋泄漏點(diǎn)更加方便快捷。泄漏參數(shù)敏感性分析發(fā)現(xiàn),管道泄漏量與泄漏位置均對(duì)管道壓力和流量信號(hào)產(chǎn)生影響。