富 友，蔣 勁 ，李燕輝，應(yīng) 銳

（1.蘭州理工大學(xué)能源與動(dòng)力工程學(xué)院，蘭州 730050；2.武漢大學(xué)水力機(jī)械過(guò)渡過(guò)程教育部重點(diǎn)實(shí)驗(yàn)室，武漢 430072）

0 引 言

當(dāng)管路中的流動(dòng)狀態(tài)發(fā)生變化時(shí)，會(huì)產(chǎn)生水錘現(xiàn)象。水錘現(xiàn)象的發(fā)生、發(fā)展程度受管路材料特性和流動(dòng)介質(zhì)屬性的影響[1]。就管路特性方面，材料的彈性模量對(duì)管路中的波速傳遞有著明顯影響，特別對(duì)于薄管壁流動(dòng)而言，水錘波的傳遞速度將會(huì)發(fā)生明顯變化[2]，對(duì)于介質(zhì)屬性而言，密度、壓力及熱力學(xué)參數(shù)均對(duì)相應(yīng)波的傳遞有明顯的影響[3]。近些年在水錘計(jì)算研究方面，主要針對(duì)水錘波的傳遞過(guò)程中受管路幾何結(jié)構(gòu)的影響，管路水錘計(jì)算的流固耦合過(guò)程及水錘壓力波的衰減過(guò)程方面進(jìn)行了較為詳細(xì)的研究[4-7]，但帶有液柱分離現(xiàn)象的水錘計(jì)算一直為計(jì)算的難點(diǎn)，有較多的學(xué)者對(duì)該種現(xiàn)象進(jìn)行過(guò)相應(yīng)的研究分析，并取得了一定的成果[8-12]。當(dāng)壓力降低至飽和蒸氣壓力以下時(shí)，由于相變導(dǎo)致的液柱分離現(xiàn)象將會(huì)對(duì)壓力響應(yīng)的準(zhǔn)確性造成更大的變化及影響[13]。

傳統(tǒng)的水錘計(jì)算方法，通常采用連續(xù)性方程和動(dòng)量方程進(jìn)行構(gòu)建計(jì)算模型，并采用特征線(xiàn)法進(jìn)行求解計(jì)算。對(duì)于計(jì)算帶有液柱分離的水錘現(xiàn)象時(shí)，通常采用離散空穴模型DVCM或DGCM模型進(jìn)行計(jì)算，計(jì)算中當(dāng)壓力降低至飽和蒸氣壓時(shí)，將所計(jì)算的區(qū)域作為計(jì)算邊界來(lái)計(jì)算相應(yīng)的空穴率的變化[14-16]。

對(duì)于DVCM模型，Simpson 等[17]采用試驗(yàn)及數(shù)值計(jì)算進(jìn)行驗(yàn)證，結(jié)果表明空穴率大于 10%時(shí)，壓力峰值及持續(xù)時(shí)間均出現(xiàn)了明顯的誤差。

Bergant等[13]指出當(dāng)空穴率大于10%時(shí)，由于交界面的問(wèn)題，將導(dǎo)致泡狀流動(dòng)下的離散空泡的交界條件不再適用。Wylie曾將一小部分自由氣體引入計(jì)算模型中，進(jìn)而阻止數(shù)值震蕩的出現(xiàn)[18]。也有學(xué)者對(duì)帶有此類(lèi)不可溶性氣體的計(jì)算方法（DGCM 模型）進(jìn)行研究，但自由氣體體積的含量具有不可確定性[1]，同時(shí)體積的量級(jí)對(duì)于計(jì)算的壓力變換有明顯的影響[19]，該方法也受到一定質(zhì)疑。不可溶性氣體的存在，使得水錘現(xiàn)象可以從兩相流的角度模擬計(jì)算[20-21]。在氣液兩相流動(dòng)計(jì)算中，雙流體模型是能夠較為準(zhǔn)確描述氣液兩相流動(dòng)特征的計(jì)算方法，該方法認(rèn)為各相的動(dòng)力學(xué)特性不完全相同，模型由 1個(gè)空泡輸運(yùn)方程、2個(gè)連續(xù)性方程、2個(gè)動(dòng)量方程及2個(gè)能量方程組成，在保證雙相的可壓縮性的前提下，方程具有良好的雙曲性[22]，可為氣液兩相瞬變流動(dòng)計(jì)算提供基礎(chǔ)。

本文選用改進(jìn)的雙流體模型進(jìn)行水錘計(jì)算模型的構(gòu)建，并通過(guò)構(gòu)建合理的相間作用形式及高精度的數(shù)值解法來(lái)描述帶有液柱分離的氣液兩相流動(dòng)現(xiàn)象。本文首先對(duì)于改進(jìn)的雙流體模型的數(shù)學(xué)特性進(jìn)行分析，并選用全空化模型對(duì)相間質(zhì)量變化進(jìn)行計(jì)算，構(gòu)建出一種考慮相間質(zhì)量變化的氣液兩相瞬變流的計(jì)算方法。其次對(duì)方程中的不守恒項(xiàng)進(jìn)行了單獨(dú)處理并選用矢通量分裂(steger-warming flux vector splitting, SW-FVS)與無(wú)振蕩，無(wú)自由參數(shù)，耗散算法 (non-oscillatory containing, no free parameters and dissipative scheme, NND)相結(jié)合的時(shí)間上為一階，空間上均為二階的數(shù)值方法進(jìn)行數(shù)值計(jì)算。最后結(jié)合2個(gè)實(shí)驗(yàn)室的試驗(yàn)數(shù)據(jù)對(duì)2種弛豫方法對(duì)于計(jì)算結(jié)果的影響進(jìn)行討論分析。

1 改進(jìn)的雙流體計(jì)算模型

水錘計(jì)算模型是基于改進(jìn)雙流體模型建立，與傳統(tǒng)的雙流體模型的區(qū)別在于方程組中添加了管壁截面參數(shù)變量，這一改變也使得該計(jì)算模型更適用于管路中的氣液兩相流動(dòng)計(jì)算。方程模型對(duì)液相和氣相分別采用不同的壓力、空穴率和密度進(jìn)行描述。模型由 1個(gè)空泡輸運(yùn)方程，2個(gè)連續(xù)性方程，2個(gè)動(dòng)量方程及2個(gè)能量方程組成。改進(jìn)的模型方程的雙曲型部分與文獻(xiàn)[3，22-23]所示相似，其中一維的數(shù)值模型建立如式（1）～式（7）。

式中αk為空穴率；ρk為密度，ρ為混合密度，ρ1為相間密度，kg/m3，k=1代表氣相，k=2代表液相；u為速度，m/s；λ為弛豫系數(shù)；Pk為壓力，Pa；E為彈性模量，MPa；e為管壁厚度，m；m ˙g為質(zhì)量傳遞率，kg/s；ck為介質(zhì)波速，ckg為管內(nèi)波速，m/s；D為直徑，m；A為面積，m2；γk為比熱比；τ為剪切應(yīng)力，N/m2；Ek為總能量，J，Ek=ek+1/2uk2；ek為內(nèi)能，J；θ 為管路傾斜角，(°)，γk為比熱比。

1.1 閉合條件

由于計(jì)算中考慮了能量的變化，因此需要引入狀態(tài)方程保證方程的封閉，這里采用修正的氣體狀態(tài)方程即SG-EOS[22]。

相間的內(nèi)能變化如式（10）～式（12）所示。

比熱比γk，參考?jí)毫Ζ衚，焓qk，定壓比熱Cvk，定容比熱Cpk參考文獻(xiàn)[23-24]。在計(jì)算中π1認(rèn)為大小為0（等同于理想氣體狀態(tài)方程），這一參數(shù)的取值帶入波速方程后可得與聲速方程一樣表達(dá)式，γk與 π2采用式（13）、（14）進(jìn)行計(jì)算。

式中Tref為溫度，K；psat為飽和蒸汽壓力，Pa。

1.2 管壁變形的計(jì)算方法

計(jì)算中僅考慮一維流動(dòng)，不考慮管道軸向的形變，僅考慮管壁徑向上形變。同時(shí)管道內(nèi)的壓力變化僅僅對(duì)管壁內(nèi)部橫截面積造成影響，對(duì)于管壁厚度不造成影響。參照文獻(xiàn)[19]，對(duì)管壁內(nèi)的橫截面積采用式（15）～（16）進(jìn)行計(jì)算。

式中 A0為初始管壁面積，m2；P∞為初始?jí)毫?，Pa，D0為初始管壁直徑，m；cn1為管壁材質(zhì)參數(shù)，v為泊松比。

2 數(shù)值計(jì)算方法

為得到更高求解精度下的計(jì)算結(jié)果，需要選用二階或更高階的求解格式進(jìn)行數(shù)值求解計(jì)算，通常對(duì)于顯示求解格式而言，高階的求解格式往往會(huì)產(chǎn)生數(shù)值震蕩，為了保證在處理壓力震蕩時(shí)的精確性及穩(wěn)定性，本文選用NND格式進(jìn)行求解計(jì)算。對(duì)于氣液兩相雙流體模型而言，方程呈現(xiàn)非守恒形式，并包含不守恒項(xiàng)，因此采用高階守恒格式的數(shù)值求解方法時(shí)，需要對(duì)相應(yīng)的不守恒方程及不守恒項(xiàng)進(jìn)行單獨(dú)處理。

2.1 空穴輸運(yùn)方程的數(shù)值計(jì)算方法

方程(1)而言，由于無(wú)法將其寫(xiě)成守恒型格式，因此需要對(duì)空泡的輸運(yùn)方程進(jìn)行單獨(dú)處理，參照修正的Godunov scheme方法[22]，對(duì)方程（1）進(jìn)行離散，其離散形式表現(xiàn)為式（17）～（18）。

2.2 雙曲型算子計(jì)算方法

對(duì)于方程（1）～（7）而言，方程中守恒部分采用SW-FVS與NND相結(jié)合的方法進(jìn)行計(jì)算。

NND是張涵信院士提出的一種數(shù)值計(jì)算方法[25]。其計(jì)算公式如式（21）～（24）所示。

2.3 對(duì)于不守恒項(xiàng)的數(shù)值計(jì)算方法

將方程（2）～（3）、（4）～（5）依照式（21）寫(xiě)成相應(yīng)的離散格式。

將方程（25）帶入方程（26）中，根據(jù)文獻(xiàn)[26]，可得表達(dá)式（28）與（29）。此時(shí)式（29）即為方程的不守恒項(xiàng)提供了相應(yīng)的求解方法。

2.4 相間弛豫的計(jì)算方法

對(duì)于雙流體模型下的氣液兩相瞬變流而言，氣液兩相分別具有各自的狀態(tài)方程，并以此來(lái)描述各相的流動(dòng)特性的變化。當(dāng)產(chǎn)生瞬變波動(dòng)時(shí)，從波速方程中可以看出，不同相間將會(huì)產(chǎn)生不同量級(jí)的波[27]。

隨著壓力波的傳遞，兩相必將會(huì)產(chǎn)生相異的壓力。對(duì)于流動(dòng)而言，介質(zhì)不可能長(zhǎng)時(shí)間保持彼此相異的環(huán)境壓力，需要一個(gè)平衡的狀態(tài)來(lái)描述該過(guò)程，同時(shí)這個(gè)過(guò)程也將導(dǎo)致兩相體積及內(nèi)能的相應(yīng)變化。通常對(duì)于兩相壓力的平衡過(guò)程是一個(gè)二維或三維的平衡過(guò)程。對(duì)于一維管道流動(dòng)而言，采用更高維度的計(jì)算方法將使計(jì)算變的更加復(fù)雜，因此需要構(gòu)建適用于一維求解的計(jì)算方法。

2.4.1 有限相間弛豫現(xiàn)象計(jì)算方法

由方程組（1）～（7）可知，方程右側(cè)存在計(jì)算壓力弛豫的源項(xiàng)，即λ（P1-P2）。對(duì)于采用兩相瞬變流動(dòng)的計(jì)算方法，可采用有限弛豫壓力系數(shù)來(lái)進(jìn)行構(gòu)建[28-30]，其構(gòu)建形式一般如式（30）形式。

式中AI為單體體積下泡粒表面積，N為泡粒數(shù)量。

這種方法需引入泡粒交界面的計(jì)算方程，同時(shí)需要考慮泡粒的脫落與融合。對(duì)于泡粒的直徑，數(shù)量密度的假設(shè)將會(huì)大大限制計(jì)算的穩(wěn)定性。并考慮到泡粒面積的改變，需要引入相間質(zhì)量的描述方式。對(duì)于交界面質(zhì)量變化的計(jì)算方法，一般認(rèn)為流動(dòng)中存在 3種物質(zhì)，即空氣、水蒸氣和水。當(dāng)壓力降低至飽和蒸氣壓力以下時(shí)，控制體內(nèi)將發(fā)生空化現(xiàn)象。選用全空化模型進(jìn)行計(jì)算相間質(zhì)量變化，全空化模型是由Singhal等[30]提出。模型考慮流動(dòng)中的不可溶性氣體，同時(shí)具有經(jīng)驗(yàn)系數(shù)少，計(jì)算準(zhǔn)確等特點(diǎn)。該種模型在CFD商業(yè)軟件中有較為廣泛的應(yīng)用，其蒸騰和凝結(jié)過(guò)程采用式 （31）～（34）計(jì)算。

式中 kt為湍動(dòng)能，m2/s2；fk為質(zhì)量分?jǐn)?shù)；Cc、Ce通常取0.01和0.02；σ為張力，N/m。

2.4.2 無(wú)限相間弛豫現(xiàn)象計(jì)算方法

當(dāng)認(rèn)為壓力弛豫現(xiàn)象的計(jì)算過(guò)程遠(yuǎn)遠(yuǎn)小于時(shí)間步長(zhǎng)時(shí)，可認(rèn)為壓力弛豫系數(shù)變換為無(wú)限大，此時(shí)可對(duì)每次計(jì)算后的結(jié)果進(jìn)行修正，以此來(lái)達(dá)到對(duì)于弛豫現(xiàn)象的描述。在使用無(wú)限弛豫現(xiàn)象修正方法的過(guò)程中，不需要采用空化模型來(lái)描述相間的質(zhì)量變化，該種方法已對(duì)相間的作用做出修正和描述。修正時(shí)源項(xiàng)對(duì)于空間的偏微分對(duì)計(jì)算不產(chǎn)生影響，參考相關(guān)文獻(xiàn)[22]有式（35）～（39）。

將式（35）帶入式（38）、（39）中，可以得到式（40）。

將式（40）進(jìn)行積分可以得到式（41）。

式中上標(biāo)0表示雙曲算子后的計(jì)算結(jié)果，上標(biāo)*表示弛豫過(guò)程后的計(jì)算結(jié)果。將方程進(jìn)行進(jìn)一步整理可得式（42）。

將狀態(tài)方程（9）、相間壓力計(jì)算方程進(jìn)行帶入，同時(shí)由于π1通常取為0，因此有式（43）。

同時(shí)從式（38）、（39）可知，過(guò)大的泡粒直徑將導(dǎo)致弛豫系數(shù)增大，進(jìn)而導(dǎo)致弛豫系數(shù)對(duì)于相間能量的反向修正[31-32]，也驗(yàn)證了2.4.1中的結(jié)論。

2.5 邊界條件

對(duì)于管道氣液兩相流動(dòng)的上、下游的邊界條件，采用特征線(xiàn)法進(jìn)行構(gòu)建，在構(gòu)建時(shí)一般需要知道上游、下游的速度或壓力的變化，以此來(lái)計(jì)算其他流動(dòng)參數(shù)的變化趨勢(shì)。其中相容性方程采用式（45）進(jìn)行插值計(jì)算。

3 數(shù)值計(jì)算及試驗(yàn)比對(duì)

為驗(yàn)證該計(jì)算模型的適用性及準(zhǔn)確性，選用了 2個(gè)實(shí)驗(yàn)室的不同試驗(yàn)數(shù)據(jù)，分別與數(shù)值模擬的結(jié)果進(jìn)行對(duì)比，試驗(yàn)包括辛普森的關(guān)閥水錘試驗(yàn)[33]、武漢大學(xué)的水錘試驗(yàn)。

3.1 辛普森關(guān)閥水錘試驗(yàn)

辛普森關(guān)閥試驗(yàn)由四部分組成，分別為上游穩(wěn)壓水池、管線(xiàn)、尾閥和下游水池組成。圖 1為辛普森關(guān)閥試驗(yàn)示意圖。

圖1 辛普森試驗(yàn)示意圖Fig. 1 Sketch of Simpson’s experiment

管線(xiàn)總長(zhǎng)為36 m，直徑為 0.019 05 m，厚度為 0.001 588 m，管道在9 m位置處有一90o轉(zhuǎn)接頭，管道末端連接低位水池，在管線(xiàn)36 m位置安裝有一球閥。管線(xiàn)具有1 m高的坡度，在管線(xiàn)沿線(xiàn)有3處壓力檢測(cè)點(diǎn)，分別位于9、27、36 m。計(jì)算初始空穴率為10-7，上游穩(wěn)壓水箱壓力為31.88 m，流速為1.125 m/s，管線(xiàn)阻力損失為2.8 m，飽和蒸氣壓力為3 032.17 Pa，波速為1 280 m/s，流體主要參數(shù)均參考文獻(xiàn)[32]。水和空氣的密度分別為998.2、2.52 kg/m3，飽和水蒸氣的密度為0.077 14 kg/m3。表1為所對(duì)應(yīng)的流體固有熱力學(xué)參數(shù)。模擬計(jì)算以尾閥處速度邊界進(jìn)行讀入，關(guān)閥時(shí)間為0.043 07 s。圖2為辛普森關(guān)閥試驗(yàn)結(jié)果與計(jì)算結(jié)果比對(duì)圖。

表1 辛普森試驗(yàn)熱力學(xué)參數(shù)Table 1 EOS parameters for liquid and gas for Simpson’s experiment

圖2 辛普森試驗(yàn)結(jié)果與計(jì)算結(jié)果比對(duì)Fig. 2 Comparison between model and experiments in Simpson’s experiment

由圖 2可知，當(dāng)采用有限弛豫系數(shù)計(jì)算方法時(shí)，如果不加入空化模型，那么計(jì)算結(jié)果與試驗(yàn)數(shù)據(jù)將出現(xiàn)極大的差異。從采用有限弛豫系數(shù)計(jì)算方法并引入空化模型的計(jì)算中可以看出，當(dāng)泡粒數(shù)量（單位體積下泡粒數(shù)）選定為 1013時(shí)，計(jì)算結(jié)果與試驗(yàn)有較好的相似度，且泡粒數(shù)量選擇量級(jí)與文獻(xiàn)[34-35]描述一致。

在計(jì)算中發(fā)現(xiàn)，泡粒數(shù)量對(duì)于計(jì)算結(jié)果有明顯的影響，當(dāng)泡粒數(shù)量小于1013時(shí)，如圖3中所示會(huì)出現(xiàn)次級(jí)波峰的響應(yīng)時(shí)間出現(xiàn)偏差。當(dāng)泡粒數(shù)量大于 1013時(shí)，會(huì)因?yàn)閱挝粷舛认碌呐萘１砻娣e過(guò)大導(dǎo)致相間的動(dòng)量傳遞超過(guò)真實(shí)計(jì)算值，造成內(nèi)能和壓力的反向修正，從而導(dǎo)致計(jì)算中斷[31]。

圖3 不同泡粒數(shù)量對(duì)辛普森試驗(yàn)的影響Fig. 3 Effect of different particle concentration in Simpson’s experiment

從圖 2計(jì)算結(jié)果中可以看出，雖無(wú)限弛豫系數(shù)的計(jì)算結(jié)果較當(dāng)泡粒數(shù)量選定為 1013時(shí)的計(jì)算結(jié)果有一定偏差，但與試驗(yàn)值的偏差同樣較小。同時(shí)通過(guò)圖2d空穴率的比對(duì)可以看出，2種計(jì)算方法計(jì)算的空穴率的變化雖然有一定偏差，但空穴率的變化趨勢(shì)相似性較高。

圖4 采用不同弛豫計(jì)算方法時(shí)氣液兩相壓力差Fig. 4 Difference of gas-liquid two-phase pressure difference with different relaxation calculation methods

圖 4為采用無(wú)限弛豫系數(shù)修正前氣相壓力與液相壓力差，與采用有限弛豫系數(shù)并選用空化模型下的壓力差對(duì)比圖。由圖 4可知，采用無(wú)限弛豫系數(shù)修正前，在壓力低于飽和蒸氣壓力的空化初生區(qū)，氣相壓力始終大于液相壓力。由式（42）可知，采用無(wú)限弛豫系數(shù)，正是對(duì)相間的能量變化進(jìn)行修正，能量的變化雖然沒(méi)有體現(xiàn)出相變，但從本質(zhì)上來(lái)說(shuō)，可以看成是伴隨著能量變化的另一種相變表達(dá)形式，因此該種修正方法下不需要再次加入相變修正模型。同樣從圖4中可以看出當(dāng)在0.3 s左右，正壓力返回至末端時(shí)，對(duì)于采用有限弛豫系數(shù)修正的方法依然造成了壓力相修正的偏差。當(dāng)選用的泡粒數(shù)量繼續(xù)增大時(shí)，將會(huì)出現(xiàn)反向修正的情況。

由于無(wú)限壓力弛豫計(jì)算方法為一種強(qiáng)制性的修正計(jì)算方法，即在每個(gè)時(shí)間步長(zhǎng)后，認(rèn)為氣相、液相壓力均達(dá)到平衡狀態(tài)，而事實(shí)上壓力的弛豫時(shí)間通常為1 ms。若計(jì)算時(shí)間步長(zhǎng)小于1 ms時(shí)，由于計(jì)算方法的強(qiáng)制修正，使得未到達(dá)平衡狀態(tài)下的兩相流動(dòng)強(qiáng)制處于平衡狀態(tài)，這也導(dǎo)致了空穴體積的修正偏大，使得在低壓區(qū)的空化模型作用被大幅度消減。

3.2 管線(xiàn)上游水錘關(guān)閥試驗(yàn)

為模擬離心泵突然失電后關(guān)閥導(dǎo)致的管路系統(tǒng)中壓力瞬變，本文在武漢大學(xué)水力機(jī)械過(guò)渡過(guò)程實(shí)驗(yàn)室，對(duì)此做了相應(yīng)試驗(yàn)進(jìn)行對(duì)比驗(yàn)證。試驗(yàn)管路系統(tǒng)布置如圖5所示。試驗(yàn)由水箱、水泵、閥門(mén)、壓力傳感器、電磁流量計(jì)及水平管線(xiàn)組成，其中試驗(yàn)段管線(xiàn)總長(zhǎng)為 97.6 m，管徑為0.1 m，管壁材料為鍍鋅鋼管，管壁厚為3 mm，折算管壁聲速約為1 260 m/s；末端為高位恒壓水箱，水位高度為 26 m。水和空氣的密度分別為 998.2和 2.52 kg/m3，飽和水蒸氣的密度為 0.077 14 kg/m3。試驗(yàn)中采取電動(dòng)蝶閥進(jìn)行快關(guān)試驗(yàn)。計(jì)算中選用無(wú)限弛豫計(jì)算方法進(jìn)行計(jì)算。試驗(yàn)數(shù)據(jù)如表2所示。

圖5 試驗(yàn)管路系統(tǒng)布置圖Fig. 5 Layout drawing of test piping system

表2 試驗(yàn)初始數(shù)據(jù)Table 2 Initial conditions of experiment

對(duì)于壓力的測(cè)試，采用閥門(mén)位置（壓力傳感器P1點(diǎn)所在位置）和0.5倍管線(xiàn)位置（壓力傳感器P2點(diǎn)所在位置）進(jìn)行數(shù)據(jù)采集。試驗(yàn)采用武漢超宇測(cè)控技術(shù)有限公司生產(chǎn)的 CY 3088型號(hào)壓力變送器，測(cè)量范圍為 0～5.0 MPa（絕對(duì)壓力），測(cè)量精度為 0.5%，流量計(jì)為KROHNE生產(chǎn)的IFM 4080K電磁流量計(jì)，測(cè)量范圍為0～80 m3/h，測(cè)量精度為0.5%，試驗(yàn)系統(tǒng)中的儀表精度等級(jí)直接決定了試驗(yàn)臺(tái)的測(cè)試精度，在本章試驗(yàn)系統(tǒng)中試驗(yàn)臺(tái)總測(cè)試系統(tǒng)的系統(tǒng)誤差為±0.933 8%。試驗(yàn)中選用的測(cè)量?jī)x器輸出均為電流信號(hào)，采用數(shù)據(jù)采集卡將電流設(shè)備轉(zhuǎn)換為數(shù)值信號(hào)，并通過(guò)傳輸存儲(chǔ)于計(jì)算機(jī)指定位置上。表 3為關(guān)閥試驗(yàn)熱力學(xué)參數(shù)。將計(jì)算結(jié)果和試驗(yàn)結(jié)果進(jìn)行比對(duì)，如圖6所示。由圖6可知，計(jì)算結(jié)果與試驗(yàn)結(jié)果吻合程度較高，驗(yàn)證了算法的可行性。

表4為試驗(yàn)對(duì)比誤差。由表4可知，試驗(yàn)1和試驗(yàn)2壓力波峰值均超過(guò)Joukowsky水錘定律（a△u/g）計(jì)算值，證明瞬變流動(dòng)中產(chǎn)生了明顯的液柱分離現(xiàn)象，同時(shí)通過(guò)壓力峰值與試驗(yàn)值的計(jì)算偏差可以看出計(jì)算值與試驗(yàn)值偏差較小，驗(yàn)證了數(shù)值方法的準(zhǔn)確性。

表3 關(guān)閥試驗(yàn)熱力學(xué)參數(shù)Table 3 EOS parameters for liquid and vapor water for experiment

圖6 關(guān)閥試驗(yàn)結(jié)果與計(jì)算結(jié)果對(duì)比Fig. 6 Comparison between experimental data and calculation data in water hammer experiment

表4 試驗(yàn)對(duì)比誤差Table 4 Error comparison between experimental and calculation

綜上 2個(gè)實(shí)驗(yàn)室的不同試驗(yàn)可以看出無(wú)限弛豫系數(shù)無(wú)空化模型的計(jì)算方法可以較好的模擬管路系統(tǒng)中，由于關(guān)閥所造成的水力過(guò)渡的壓力瞬變，并可以較為精準(zhǔn)的預(yù)測(cè)管路系統(tǒng)的過(guò)渡過(guò)程。同時(shí)相對(duì)于以往的氣液兩相計(jì)算模型[37]而言，本文所選用的計(jì)算模型，去掉了泡狀流假設(shè)，可以使得計(jì)算的穩(wěn)定性得到更好的保障。

4 結(jié) 論

本文基于改進(jìn)的雙流體氣液兩相流計(jì)算模型，對(duì)管道中由于水錘造成的液柱分離的氣液兩相瞬變流動(dòng)進(jìn)行計(jì)算分析比對(duì)，計(jì)算中考慮了管壁截面隨壓力變化時(shí)對(duì)波速的影響，并對(duì)方程中的不守恒項(xiàng)及不守恒方程進(jìn)行了單獨(dú)處理，以確保模型計(jì)算的穩(wěn)定性。通過(guò)對(duì)比兩種弛豫計(jì)算方法對(duì)于計(jì)算的影響，并采用辛普森關(guān)閥試驗(yàn)與管線(xiàn)上游水錘關(guān)閥試驗(yàn)進(jìn)行分析比對(duì)論證，結(jié)論如下：

1）采用有限弛豫參數(shù)計(jì)算時(shí)，需考慮由于空化引起的相變，同時(shí)需選擇合理的泡粒數(shù)量以保證計(jì)算精確度。

2）采用無(wú)限弛豫參數(shù)計(jì)算時(shí)，空化模型對(duì)于計(jì)算結(jié)果影響較小，可直接對(duì)水錘現(xiàn)象造成的液柱分離現(xiàn)象進(jìn)行分析模擬，但會(huì)造成空穴率計(jì)算結(jié)果過(guò)大。

3）采用改進(jìn)的雙流體氣液兩相流計(jì)算模型，并構(gòu)建合理的數(shù)值計(jì)算方法及弛豫系數(shù)，可以對(duì)帶有明顯液柱分離現(xiàn)象的管路氣液兩相瞬變現(xiàn)象進(jìn)行預(yù)測(cè)計(jì)算，該方法計(jì)算精度高，對(duì)于工業(yè)管道中的水錘防護(hù)計(jì)算具有較好的應(yīng)用價(jià)值。