彈性/黏彈性混合管道瞬變流參數(shù)校核方法

孫強(qiáng)， 王祉皓， 伍悅濱， 徐瑩

(1.東北林業(yè)大學(xué)土木工程學(xué)院人工環(huán)境與能源應(yīng)用研究所， 哈爾濱 150040; 2.哈爾濱工業(yè)大學(xué)建筑學(xué)院， 哈爾濱 150090； 3.哈爾濱工業(yè)大學(xué)寒地建筑科學(xué)與工程研究中心， 哈爾濱 150090; 4.哈爾濱商業(yè)大學(xué)能源與建筑工程學(xué)院， 哈爾濱 150028)

有壓管道經(jīng)常會因突然停泵或關(guān)閥引起的壓力突增或驟降的水力現(xiàn)象，可能導(dǎo)致設(shè)備損壞或管道炸裂，這種管內(nèi)流體的瞬變流動現(xiàn)象，在工程中稱之為水錘現(xiàn)象[1-2]。為保證輸水系統(tǒng)安全運(yùn)行，對管道瞬變流動進(jìn)行分析和控制是學(xué)者們重點(diǎn)研究的方向[3-4]。

近些年，隨著高分子材料技術(shù)的迅猛發(fā)展，高分子聚合物管道因其價格低、抗腐蝕性高、節(jié)能環(huán)保等特點(diǎn)，廣泛應(yīng)用于供水管道系統(tǒng)中?，F(xiàn)如今，塑料管道常用于替換傳統(tǒng)金屬管道，使得現(xiàn)階段管道系統(tǒng)中經(jīng)常出現(xiàn)金屬管道與塑料管道混合連接的情況。由于塑料管道的力學(xué)特性不同于傳統(tǒng)金屬管道的彈性力學(xué)特性，其同時表現(xiàn)出彈性和黏性兩種力學(xué)行為，因而也被稱為黏彈性管道。黏彈性管道在瞬變應(yīng)力的激勵下，所產(chǎn)生的應(yīng)變可以表示為彈性應(yīng)變和延遲應(yīng)變的疊加。這種黏彈力學(xué)行為還會影響瞬變流動中壓力波動的峰值、相位和衰減。因黏彈性管道力學(xué)行為的復(fù)雜性，傳統(tǒng)彈性管道瞬變流的數(shù)學(xué)模型已不再適用于黏彈性管道瞬變流的數(shù)值模擬。近些年來，學(xué)者開始對黏彈性管道瞬變流壓力波動的影響機(jī)制進(jìn)行了理論和實(shí)驗研究[5-6]，并先后提出并優(yōu)化了黏彈性管道瞬變流計算模型[7-9]，證明了這些模型的準(zhǔn)確性。

在黏彈性管道瞬變流模型中，管道本構(gòu)模型主要采用廣義開爾文-福伊特(Kelvin-Voigt，KV)力學(xué)模型，該模型中包含蠕變?nèi)崃亢脱舆t時間兩個參數(shù)，蠕變參數(shù)對瞬變流壓力波動的影響規(guī)律分析和參數(shù)在本構(gòu)模型中的選取一直是學(xué)者們重點(diǎn)關(guān)注的問題所在。Urbanowicz等[10]分析了不同水溫下蠕變參數(shù)對黏彈性管道瞬變流壓力波曲線的影響規(guī)律。Javadi等[11]分析了黏彈性管道瞬變流壓力波動周期對蠕變參數(shù)的敏感性。為探究蠕變參數(shù)的校核方法，Covas等[5]利用力學(xué)拉伸實(shí)驗，測得高密度聚乙烯(high denisity polyethylene, HDPE)管道的本構(gòu)曲線，并將其用于數(shù)值模擬，結(jié)果表明模擬結(jié)果與實(shí)驗結(jié)果偏差較大。Pezzinga[12]和Pan等[13-14]分別應(yīng)用遺傳算法和頻域分析方法對黏彈性管道瞬變流模型中的蠕變參數(shù)進(jìn)行校核。但這些參數(shù)校核方法存在校核程序復(fù)雜，耗時長等問題。

對于管道瞬變流動的研究，前人主要針對單一管材的管道系統(tǒng)，而針對彈性/黏彈性混合管道系統(tǒng)的研究較少。為闡明混合管道的瞬變流壓力波動影響規(guī)律，Garg等[15]進(jìn)行了鋼管/玻璃纖維塑料混合管道的瞬變流實(shí)驗，在實(shí)驗中改變塑料管段長度在整個管道系統(tǒng)中的占比，從而分析混合管道瞬變流壓力波動影響規(guī)律。Trabelsi等[16]在彈性管道系統(tǒng)上游設(shè)置塑料短管，研究塑料短管對水錘抑制作用。但目前還未提出適用于彈性/黏彈性混合管道的高效且準(zhǔn)確的參數(shù)校核方法。

因此，現(xiàn)搭建彈性/黏彈性混合管道瞬變流實(shí)驗臺，并進(jìn)行了鋼管、鋼管/PPR混合管道、鋼管/HDPE混合管道的瞬變流實(shí)驗，分析混合管道瞬變流壓力波動的影響機(jī)制，闡明了波速和蠕變參數(shù)對壓力波動的影響規(guī)律，進(jìn)而提出一種簡單且準(zhǔn)確的混合管道波速和蠕變參數(shù)的校核方法，基于實(shí)驗數(shù)據(jù)驗證該方法的準(zhǔn)確性，為工程中混合管道的模擬計算和參數(shù)校核提供參考和依據(jù)。

1 數(shù)學(xué)模型及求解方法

1.1 數(shù)學(xué)模型

經(jīng)典的彈性管道瞬變流控制方程包括連續(xù)性方程和動量方程[17]，分別表示為

(1)

(2)

式中：H為壓力水頭；V為平均流速；a為波速；g為重力加速度；x為距離；t為時間；JQ為擬穩(wěn)態(tài)摩阻引起的水頭損失，可以使用達(dá)西-魏斯巴赫(Darcy-Weisbach)公式為

(3)

式(3)中：Q為流量；f為沿程阻力系數(shù)；D為管徑；A為管道斷面截面積。

管道瞬變流壓力波的波速a計算公式為

(4)

式(4)中：E為楊氏模量；KL為水的體積模量；ρ為水的密度；e為管道壁厚；C0為約束系數(shù)，其計算公式為

(5)

式(5)中：μ為管道泊松比。

相比彈性管道，黏彈性管道瞬變流控制方程中的連續(xù)性方程中增加了延遲應(yīng)變項，用以模擬管道黏彈性效應(yīng)，其連續(xù)性方程為

(6)

式(6)中：εr為管道的延遲應(yīng)變。

黏彈性管道的總應(yīng)變ε為

ε=εe+εr

(7)

式(7)中：εe為黏彈性材料的瞬時應(yīng)變。

根據(jù)玻爾茲曼疊加原理，總應(yīng)變ε可表示為

(8)

式(8)中：J0為瞬時蠕變?nèi)崃浚籎(t′)為t′時刻的蠕變?nèi)崃?；?t)為t時刻的應(yīng)力。

黏彈性材料的蠕變函數(shù)使用廣義KV模型表式為

(9)

式(9)中：J0為瞬時蠕變?nèi)崃?；JK為第K個KV元件的蠕變?nèi)崃?；τk為延遲時間；n為蠕變元件個數(shù)。

1.2 模型求解

(1)彈性管道瞬變流特征線法求解。使用特征線法，將連續(xù)性方程[式(1)]和動量方程[式(2)]轉(zhuǎn)化為常微分方程，稱之為相容性方程，即

(10)

將式(10)分別沿正負(fù)特征線(圖1中C+和C-特征線)對相容性方程進(jìn)行積分，得到代數(shù)方程為

圖1 差分網(wǎng)格Fig.1 Grid of characteristics

(11)

(2)黏彈性管道瞬變流特征線法求解。使用特征線法，將黏彈性管道連續(xù)性方程[式(6)]和動量方程[式(2)]轉(zhuǎn)化為常微分方程，即

(12)

分別沿正負(fù)特征線對方程進(jìn)行積分，得到代數(shù)方程為

(13)

2 實(shí)驗裝置與結(jié)果分析

2.1 實(shí)驗設(shè)備及方法

如圖2所示，彈性/黏彈性混合管道瞬變流實(shí)驗臺是由供水水箱、管道、氣動閥、替換活接、下游水槽等組成。實(shí)驗管段為總長為15 m的鍍鋅鋼管，管道的公稱直徑為DN25、壁厚為3 mm，彈性模量為1.4×1014Pa，泊松比μ為0.23。上游水箱的水位為2.7 m。在管道末端設(shè)置氣動球閥作為管道瞬變流的激勵方式。管道中段兩處活接的間距為5.7 m，活接部分管道占整個管道的38%，關(guān)閥時間0.15 s，線性關(guān)閥，通過實(shí)驗數(shù)據(jù)分析和理論計算確定管道沿程阻力系數(shù)f為0.036。實(shí)驗設(shè)備的具體型號參數(shù)，如表1所示。實(shí)驗裝置的局部實(shí)物圖如圖3所示。

圖2 實(shí)驗裝置示意圖Fig.2 Schematic diagram of experimental setup

表1 設(shè)備及主要參數(shù)Table 1 Main parameters of equipment

本次實(shí)驗分兩部分進(jìn)行，第一部分為彈性管道瞬變流實(shí)驗，通過關(guān)閉管道下游的氣動閥門產(chǎn)生管道瞬變流動；第二部分為彈性與黏彈性混合管道瞬變流實(shí)驗，使用管道中預(yù)先設(shè)置好的活接，將管道中段鋼管替換為高密度聚乙烯(HDPE)和聚丙烯(PPR)管道，從而進(jìn)行彈性/黏彈性混合管道瞬變流實(shí)驗，所替換的黏彈性管段的管徑和壁厚與彈性管道相同。所進(jìn)行的實(shí)驗工況，如表2所示。

本次實(shí)驗所采用的壓力傳感器的精度為±0.25% FS，在本實(shí)驗中測得壓力波的峰值為12 m，滿度相對誤差為1.4%，在可接受的范圍內(nèi)。另外，工況2的水溫為24 ℃，其他工況的水溫均為15 ℃。

圖3 實(shí)驗裝置局部實(shí)物圖Fig.3 Partial physical picture of experimental device

表2 實(shí)驗工況基本參數(shù)Table 2 Parameters of experimental tests

2.2 實(shí)驗數(shù)據(jù)分析

對上游水箱作為邊界條件、末端關(guān)閥作為瞬變流激勵方式的實(shí)驗進(jìn)行分析。在相同實(shí)驗參數(shù)(表2)下，分別進(jìn)行鍍鋅鋼管道(工況1)、鋼管/PPR管道(工況4)和鋼管/HDPE(工況5)管道瞬變流壓力波動實(shí)驗數(shù)據(jù)，如圖4所示。

從圖4可見，彈性管道瞬變流壓力波的最大波峰值為11.3 m，而鋼管/PPR混合管道瞬變流壓力波的最大波峰值為10.2 m，鋼管/HDPE混合管道瞬變流壓力波最大波峰值為8.8 m。與彈性管道瞬變流壓力波曲線相比，鋼管/PPR混合管道和鋼管/HDPE混合管道的壓力波最大峰值分別下降9%和22%。另外，在部分管段替換為HDPE管材后，其產(chǎn)生的瞬變流壓力波曲線相比彈性管道會出現(xiàn)更加嚴(yán)重的衰減現(xiàn)象，但在替換為PPR管道后衰減幅度與彈性管道相似。但是，鋼管/PPR混合管道相比鋼管/HDPE混合管道產(chǎn)生了更大的相位延遲現(xiàn)象。

圖4 工況1、工況4和工況5的實(shí)驗數(shù)據(jù)比較Fig.4 Comparison of experimental data for Cases 1， 4 and 5

3 參數(shù)校核方法及驗證

提出了一種合理且簡便的彈性/黏彈性混合管道瞬變流參數(shù)校核方法，并基于不同工況下的實(shí)驗數(shù)據(jù)，利用黏彈性管道瞬變流模型驗證該校核方法的準(zhǔn)確性。對于黏彈性管道，非穩(wěn)態(tài)摩阻對壓力波動的影響較小，而黏彈性效應(yīng)對管道瞬變流壓力波動的影響占主導(dǎo)地位[18-19]，因此本文研究中將使用黏彈性管道瞬變流模型結(jié)合一維擬穩(wěn)態(tài)摩阻模型進(jìn)行校核。

3.1 彈性管道瞬變流波速校核

以工況1的實(shí)驗數(shù)據(jù)為例，對彈性管道瞬變流波速進(jìn)行校核。通過實(shí)驗數(shù)據(jù)計算波速的范圍區(qū)間，基于實(shí)驗數(shù)據(jù)計算波速的方法有以下兩種。

(1)在下游氣動閥瞬時關(guān)閉后，壓力波依次通過壓力傳感器2和壓力傳感器1，分別記錄兩傳感器首次出現(xiàn)壓升的時間、壓力首次達(dá)到波峰15%的時間、壓力首次達(dá)到波峰的時間。根據(jù)公式a=(L1-L2)/(T1-T2)，可計算出壓力波波速。

(2)基于壓力波每次到達(dá)穩(wěn)態(tài)的時間計算出壓力波的周期值，根據(jù)公式a=4L/T，計算出壓力波波速。

應(yīng)用所述兩種方法計算出工況1的波速均為350 m/s左右，考慮到實(shí)驗數(shù)據(jù)存在誤差，仍需要應(yīng)用管道瞬變流模型對波速進(jìn)行進(jìn)一步校核。

根據(jù)式(4)計算彈性管道的瞬變流壓力波波速，水體積模量KL取2.1×109Pa、空氣體積模量Kg取 1×105Pa、鋼管彈性模量E取1.4×1011Pa。泊松比取0.23。根據(jù)式(5)計算出約束系數(shù)C0為1.17，計算瞬變流波速為1 100 m/s。理論計算結(jié)果同實(shí)驗結(jié)果差距過大，推測為管道含氣的原因。根據(jù)文獻(xiàn)[20]中的含氣管道波速計算公式，當(dāng)波速為350 m/s時，管道含氣量為0.000 74。

通過實(shí)驗數(shù)據(jù)計算出波速為350 m/s后，以實(shí)驗波速的±5%作為取值范圍，通過彈性管道瞬變流模型進(jìn)一步校核波速。在彈性管道瞬變流動過程中，波速是影響壓力波曲線變化的決定性因素，直接影響壓力波的周期和首個波峰值。在取值范圍內(nèi)改變波速，找到壓力波首個波峰和周期與實(shí)驗曲線契合程度最好的波速值即為校核結(jié)果。波速值為341 m/s為工況1的校核結(jié)果，模擬結(jié)果如圖5所示。結(jié)果表明，使用一維擬穩(wěn)態(tài)摩阻模型可以較為準(zhǔn)確的模擬出壓力波動的峰值和周期情況，但無法準(zhǔn)確表達(dá)壓力波動的衰減情況，這是由于彈性管道的瞬變流壓力波動的衰減很大程度上受到非穩(wěn)態(tài)摩阻的影響，而一維擬穩(wěn)態(tài)摩阻模型低估了瞬變流壓力波動過程的摩阻耗散[21]。

圖5 工況1的實(shí)驗數(shù)據(jù)及模擬結(jié)果對比Fig.5 Comparison of experimental and simulation results in case 1

3.2 混合管道水錘波速校核

用3.1節(jié)相同的方法對工況2的實(shí)驗數(shù)據(jù)進(jìn)行計算，計算得出工況2的實(shí)驗波速的估算值為230 m/s。與彈性管道相同，實(shí)驗計算波速并非準(zhǔn)確的結(jié)果，需要應(yīng)用黏彈性管道瞬變流模型進(jìn)行進(jìn)一步校核。

在參數(shù)校核的過程中，瞬變流波速大小的影響主要體現(xiàn)在壓力波的峰值和周期。但壓力波的峰值和周期同時受到波速和蠕變參數(shù)的影響，僅僅使用波峰作為波速的校核基準(zhǔn)無法得到準(zhǔn)確的結(jié)果。在瞬變流壓力波動過程中前3個周期，管壁的黏彈性對壓力波周期的影響較小，可以忽略不計，但在第4個壓力波動周期以后，受管壁的黏彈性影響，壓力波會產(chǎn)生少量的相位延遲現(xiàn)象。所以，使用壓力波動的前3個周期作為基準(zhǔn)可以較準(zhǔn)確校核出壓力波速值。以本實(shí)驗為例，在黏彈性管道瞬變流模型中選取波速值為251 m/s時，模擬結(jié)果的前3個周期與實(shí)驗曲線的契合程度最高。

3.3 本構(gòu)參數(shù)校核

對于黏彈性管道的數(shù)值模擬，僅僅通過波速不能充分模擬出壓力波動情況，還需要在本構(gòu)模型(即KV模型)中確定蠕變參數(shù)。并且，由于混合管道中管材占比不同，通過力學(xué)實(shí)驗直接測量管道的蠕變參數(shù)的方法是不可行的，只能借助實(shí)驗數(shù)據(jù)使用參數(shù)校核的方式得到混合管道的蠕變參數(shù)。

黏彈性管道本構(gòu)模型(即KV模型)常用于表達(dá)黏彈性管道的蠕變動力學(xué)行為。廣義的KV模型是多個KV元件組成的，每個KV元件是由一個彈簧和粘壺并聯(lián)組成的。KV元件的個數(shù)影響蠕變曲線的準(zhǔn)確性，選擇2個及以上的KV元件即可準(zhǔn)確模擬出管道的黏彈性行為[22]。蠕變參數(shù)的影響體現(xiàn)在壓力波的衰減和延遲上。在模型中，蠕變參數(shù)由蠕變?nèi)崃亢脱舆t時間表達(dá)。在參數(shù)校核時將延遲時間固定在0.05 s和0.25 s，通過改變?nèi)渥內(nèi)崃渴沟媚M曲線接近實(shí)驗曲線，確定出使模擬結(jié)果最接近實(shí)驗曲線的蠕變參數(shù)。

表3列出3組不同大小的蠕變參數(shù)，蠕變參數(shù)的大小對壓力波曲線的影響如圖6所示，其中a組的蠕變參數(shù)為工況2校核結(jié)果，模擬曲線同實(shí)驗曲線契合程度最高。使用該組蠕變參數(shù)以及3.2節(jié)校核得出得瞬變流波速進(jìn)行模擬計算，模擬結(jié)果如圖7所示。使用黏彈性管道瞬變流模型結(jié)合一維擬穩(wěn)態(tài)摩阻模型可以較好地描述黏彈性管道瞬變流壓力波動過程，只有第2、3個周期壓力波的衰減會有所差異，這是由于在模型中沒有考慮非穩(wěn)態(tài)摩阻的影響。

表3 三組蠕變參數(shù)Table 3 Three sets of creep parameters for simulation

圖6 不同蠕變參數(shù)下PPR混合管道的模擬結(jié)果Fig.6 Simulation results of PPR mixed pipeline under different creep functions

圖7 工況2的實(shí)驗數(shù)據(jù)及模擬結(jié)果對比Fig.7 Comparison of experimental and simulation results in case 2

如圖6所示，在延遲時間固定的情況下，b組中使用更大的蠕變參數(shù)會導(dǎo)致更大的壓力波衰減，使得壓力波的過早得趨于平緩，同時相位也會產(chǎn)生更大的延遲現(xiàn)象，c組同理。在校核過程中，第一個波峰值也會在一定程度上受到蠕變參數(shù)的影響，尤其在關(guān)閥時間較長的系統(tǒng)中，蠕變參數(shù)對第一個峰值的影響更為明顯。所以，在校核蠕變參數(shù)的過程中應(yīng)綜合考慮波峰值、衰減幅度和相位延遲情況。

3.4 不同流量的鋼管/PPR混合管道校核模擬

以工況3和工況4的實(shí)驗數(shù)據(jù)為基礎(chǔ)，對不同流量的工況進(jìn)行模擬校核并驗證不同流量下數(shù)學(xué)模型的準(zhǔn)確性。該兩組實(shí)驗在同溫度、同約束條件下進(jìn)行，工況3、工況4的初始流量分別為0.48、0.58 m3/h。

利用前文所提到的主要參數(shù)校核方法，校核得出該兩組實(shí)驗的波速為303 m/s，延遲時間τ1和τ2分別為0.05、0.25 s，蠕變?nèi)崃縅1和J2分別為0.19×10-9、0.6×10-9Pa-1。工況3和工況4的波速與本構(gòu)參數(shù)同工況2下稍有不同的原因是因為受到了水溫等外界條件影響。更高的水溫會導(dǎo)致更大的蠕變參數(shù)和更小的波速。

使用校核得出的主要參數(shù)進(jìn)行模擬的結(jié)果，如圖8所示。從圖8可知，使用黏彈性管道瞬變流模型結(jié)合一維擬穩(wěn)態(tài)摩阻模型可以在一定程度上較準(zhǔn)確地模擬出不同初始流量下的混合管道瞬變流壓力波動，模擬結(jié)果與實(shí)驗結(jié)果產(chǎn)生一定相位偏差是在可接受的范圍之內(nèi)。

圖8 工況3和工況4的實(shí)驗數(shù)據(jù)及模擬結(jié)果對比Fig.8 Comparison of experimental and simulation results in case 3 and case 4

3.5 鋼管/HDPE混合管道的校核模擬

基于工況5的鋼管/HDPE混合管道瞬變流實(shí)驗數(shù)據(jù)，利用所提到的參數(shù)校核方法進(jìn)行數(shù)值模擬。校核得出該工況下的波速為290 m/s、延遲時間τ1和τ2分別為0.05 s和0.5 s，蠕變?nèi)崃縅1和J2分別為0.13×10-9Pa-1和0.5×10-9Pa-1。模擬結(jié)果與實(shí)驗結(jié)果的對比，如圖9所示。使用黏彈性管道瞬變流模型結(jié)合一維擬穩(wěn)態(tài)摩阻模型可以較準(zhǔn)確地模擬出混合管道瞬變流壓力波動前期的峰值和衰減情況，誤差可以保持在4%以下。但模擬結(jié)果與實(shí)驗結(jié)果的相位出現(xiàn)一定程度的偏差，達(dá)到10%，在可以接受的范圍之內(nèi)。

4 結(jié)論

通過設(shè)計并搭建彈性/黏彈性混合管道瞬變流實(shí)驗臺，分別進(jìn)行了鍍鋅鋼管、鋼管/PPR混合管道以及鋼管/HDPE混合管道的瞬變流實(shí)驗，分析彈性/黏彈性混合管道的瞬變流壓力波動機(jī)制。同時提出了一種簡單準(zhǔn)確的波速和蠕變參數(shù)校核方法，并應(yīng)用黏彈性管道瞬變流模型結(jié)合一維擬穩(wěn)態(tài)摩阻模型對該方法的準(zhǔn)確性進(jìn)行驗證。主要得到如下結(jié)論。

圖9 工況5的實(shí)驗數(shù)據(jù)及模擬結(jié)果對比Fig.9 Comparison of experimental and simulation results in case 5

(1)在彈性管道系統(tǒng)中部分替換為黏彈性管道可以使得瞬變流壓力波波速變小、峰值下降及相位延遲。這一現(xiàn)象在鋼管/PPR混合管道中相比鋼管/HDPE混合管道更為明顯，在本文的實(shí)驗中，對于黏彈性管道占比為38%的混合管道，鋼管/PPR混合管道和鋼管/HDPE混合管道的壓力波最大峰值分別下降9%和22%。

(2)使用彈性管道一維擬穩(wěn)態(tài)摩阻模型無法準(zhǔn)確地模擬出彈性管道瞬變流壓力波動衰減情況。但是不同于彈性管道，使用黏彈性管道瞬變流模型結(jié)合一維擬穩(wěn)態(tài)摩阻模型可以較準(zhǔn)確地對混合管道瞬變流壓力波動進(jìn)行的數(shù)值模擬，盡管在衰減和相位上有一定的偏差，但均在可接受的范圍之內(nèi)。

(3)提出了黏彈性管道瞬變流壓力波波速和蠕變參數(shù)的校核方法。該方法首先通過實(shí)驗數(shù)據(jù)計算波速的估算區(qū)間。其次，在此區(qū)間內(nèi)以壓力波前3個周期為基準(zhǔn)，校核出與實(shí)驗曲線的前3個周期契合程度最高的波速值。最后，以壓力波的首個波峰和周期作為基準(zhǔn)，校核出與實(shí)驗曲線的峰值和相位契合度最高的蠕變參數(shù)值。