肖鵬飛，趙志高，楊桀彬，劉程鵬，何相慧，楊建東

［1．中國(guó)長(zhǎng)江三峽集團(tuán)有限公司科學(xué)技術(shù)研究院，北京市 101100；2．水資源工程與調(diào)度全國(guó)重點(diǎn)實(shí)驗(yàn)室（武漢大學(xué)），湖北省武漢市 430072］

0 引言

隨著風(fēng)光新能源的發(fā)展，抽水蓄能成為新型電力系統(tǒng)中的基礎(chǔ)性調(diào)節(jié)電源。根據(jù)《抽水蓄能中長(zhǎng)期發(fā)展綱要》，中國(guó)抽水蓄能裝機(jī)將翻兩番，達(dá)到1.2 億kW。在資源儲(chǔ)備站點(diǎn)日益緊張的背景下，抽水蓄能向著高水頭、大容量、長(zhǎng)輸水管道方向發(fā)展。抽水蓄能機(jī)組復(fù)雜的流道布置、頻繁的工況切換[1]可能導(dǎo)致水力系統(tǒng)出現(xiàn)較大的水錘振動(dòng)、液柱分離甚至空化空蝕現(xiàn)象[2]。如何準(zhǔn)確模擬尾水管液柱分離對(duì)判別抽水蓄能運(yùn)行安全性、確定機(jī)組安裝高程具有重要意義。

國(guó)內(nèi)外研究學(xué)者在20 世紀(jì)80 年代就對(duì)伴有空穴流、液柱分離的氣液兩相瞬態(tài)過程展開了研究。楊建東[3]從空化和氣液兩相流理論兩個(gè)方面，采用理論與實(shí)驗(yàn)相結(jié)合的方法，系統(tǒng)地研究了空化流動(dòng)與液核分離的瞬態(tài)過程，為今后進(jìn)一步研究其瞬態(tài)過程的機(jī)理、數(shù)值分析和滿足實(shí)際需要奠定了基礎(chǔ)。蔣勁等[4]給出了一種考慮兩相間物質(zhì)質(zhì)量交換、管壁流質(zhì)相互作用的一維計(jì)算模型。周領(lǐng)等[5]將Wylie[6]提出的離散空腔模型（DGCM）和有限體積法Godunov 格式結(jié)合起來計(jì)算水柱分離。楊桀彬[7]在集中空穴模型和氣泡均勻分布模型的基礎(chǔ)上提出了以臨界空穴率作為水柱分離的充分條件的氣液兩相流計(jì)算模型。何相慧[8]在抽水蓄能機(jī)組動(dòng)態(tài)特性實(shí)驗(yàn)平臺(tái)上開展了水柱分離模型實(shí)驗(yàn)，給出了尾水管發(fā)生水柱分離的判別條件。

本文基于有限體積法（FVM），引入空穴率α描述氣體體積，采用二階Godunov 格式離散雙曲型偏微分方程，對(duì)含氣型氣液兩相流以及相變下的瞬態(tài)過程進(jìn)行了數(shù)值模擬，尤其是摒棄了Wylie 離散空腔模型，直接采用空化模型封閉方程組，并通過工程實(shí)例與傳統(tǒng)計(jì)算方法進(jìn)行對(duì)比。

1 氣液兩相流及空化模型

1.1 基于有限體積法的氣液兩相流模型

為了衡量管道中含氣量變化，引入量度氣液兩相流中氣泡含量大小的參數(shù)，容積空穴率α，其表達(dá)式為：

式中：W——體積；

g和l——?dú)庀嗪鸵合啵?/p>

L——管道長(zhǎng)度；

A——段面橫截面面積。

氣液兩相瞬變流的基本方程仍然是動(dòng)量方程和能量方程，應(yīng)對(duì)氣相和液相分別列出。對(duì)于含氣型氣液兩相瞬變流中氣相體積較?。é?0.15）的情況，氣泡會(huì)跟隨液體一起運(yùn)動(dòng)，于是忽略動(dòng)量交換，將動(dòng)量方程合并為一個(gè)。而連續(xù)方程仍保留分相的形式，以便計(jì)算氣液兩相之間的質(zhì)量交換[9]。

類似單相流體動(dòng)量方程推導(dǎo)方法，可得到含氣型氣液兩相瞬變流動(dòng)量方程：

式中：p——管道中心點(diǎn)壓強(qiáng)，Pa；

V——流速，m/s；

x——從任意點(diǎn)開始的沿管道軸線的坐標(biāo)距離，m；

t——時(shí)間，s；

α——空穴率；

ρg——?dú)庀嗝芏龋琸g/m3；

ρl——液相密度，kg/m3；

D——管道直徑，m；

f——達(dá)西-威斯巴哈（Darcy-Weisbach）系數(shù)；

g——重力加速度，m/s2。

連續(xù)性方程的本質(zhì)是質(zhì)量守恒方程，為考慮兩相之間的質(zhì)量交換，將兩相間的質(zhì)量變換率=dm/dt引入單相流的連續(xù)方程，可以分別寫出氣相和液相的連續(xù)方程：

對(duì)上述方程進(jìn)行變形再通過有限體積法進(jìn)行離散，得到計(jì)算表達(dá)式：

式中：Fi+12(UL,i+12,UR,i+12)，F(xiàn)i-12(UL,i-12,UR,i-12)——界面i+1/2，i-1/2的數(shù)值通量；

UL和UR——界面左右的變量值。

1.2 變波速模型與相變模型

在參考文獻(xiàn)[9]中，給出了波速方程：

式中：a——水錘波波速，m/s；

C1——管道固定方式系數(shù)；

E——管道彈性模擬量；Kg和Kl——?dú)怏w容變彈性模擬量和液體容變彈性模量。

對(duì)于空化現(xiàn)象的數(shù)值模擬，需要定量描述氣液之間的物理相變過程，關(guān)鍵在于建立或選取合適的空化模型。空化模型是用于描述氣相和液相之間質(zhì)量交換的數(shù)學(xué)模型，目前應(yīng)用較為廣泛的空化模型主要分為兩類：

（1）基于正壓流體狀態(tài)方程的空化模型；

（2）基于質(zhì)量輸運(yùn)的空化方程[10]。

本文采用的計(jì)算模型需要定量考慮氣液兩相之間的相變關(guān)系，故選取基于質(zhì)量輸運(yùn)的空化模型。有壓管道的空化現(xiàn)象，往往伴隨著氣泡的破碎和融合，在極短的時(shí)間內(nèi)就會(huì)發(fā)展為空化空腔，因而基于單個(gè)氣泡的空化模型也不適用于方程的計(jì)算和求解。最后基于相變物理本質(zhì)，選取了以Singhal 模型為基礎(chǔ)的空化模型進(jìn)行數(shù)值模擬分析，并根據(jù)氣液兩相管流實(shí)際情況對(duì)公式參數(shù)進(jìn)行了優(yōu)化。

式中：Vch——特征速度，它反映了液體和蒸汽之間局部相對(duì)速度的影響，通常認(rèn)為大空隙不移動(dòng)，所以Vch=V；

σ——表面張力系數(shù)。

在該模型中，特征速度、表面張力系數(shù)、不可凝結(jié)氣體的影響也被考慮在內(nèi)，因而被稱為 “全空化模型”。

2 有限體積法和特征線法的耦合

2.1 計(jì)算模型

在抽水蓄能管道系統(tǒng)中，可能發(fā)生空化的區(qū)域較少，實(shí)際上一般只發(fā)生在尾水管進(jìn)口。其中尾水管發(fā)生的空化和液柱分離會(huì)對(duì)水輪機(jī)和水力系統(tǒng)造成極大影響。在實(shí)際工程中，往往降低安裝高程以防止尾水管內(nèi)的壓力降至蒸汽壓力以下[11]。而需要重點(diǎn)關(guān)注的也是尾水管部分的空化現(xiàn)象。

因此，同時(shí)考慮尾水管部分空化模擬的準(zhǔn)確性和計(jì)算效率，提出了兩個(gè)解決方案。

（1）將FVM 和MOC 方法進(jìn)行耦合，采用FVM 計(jì)算尾水管，MOC 方法計(jì)算抽水蓄能電站其余水力單元。

（2）尾水管的FVM 模擬計(jì)算和其他部分的MOC 模擬計(jì)算采用不同的時(shí)間步長(zhǎng)進(jìn)行耦合。

基于Riemann 不變量，有限體積法在邊界上有以下方程：

結(jié)合特征線法基本邊界方程：

如圖1 所示，對(duì)于MOC 下游邊界和FVM 上游邊界，也就是邊界1 耦合，通過式（9）、式（10）以及邊界方程聯(lián)立可求解。類似地，對(duì)于FVM 下游邊界和MOC 上游邊界，也就是邊界2 耦合，通過式（8）、式（11）以及邊界方程聯(lián)立可求解。

圖1 邊界耦合示意圖Figure 1 Schematic diagram of boundary coupling

2.2 方法驗(yàn)證

如圖2 所示，采用串聯(lián)點(diǎn)作為邊界耦合點(diǎn)，對(duì)于上庫-閥門簡(jiǎn)單的管道系統(tǒng)，壓力管道L1 和壓力管道L2 的長(zhǎng)度均為1000m，糙率n=0.014。

圖2 水庫-閥門系統(tǒng)Figure 2 Reservoir-valve system

上游庫水位20m，閥門流量由4m3/s 在0 時(shí)刻減為0m3/s。

對(duì)于管道L1 采用FVM 進(jìn)行數(shù)值模擬計(jì)算，對(duì)于管道L2采用MOC 進(jìn)行數(shù)值模擬計(jì)算。方案一FVM 和MOC 時(shí)間步長(zhǎng)均為0.04s；方案二FVM 時(shí)間步長(zhǎng)為0.04s，MOC 時(shí)間步長(zhǎng)為0.2s。

閥門處測(cè)壓管水頭隨時(shí)間變化值如圖3 所示。顯然，在不同時(shí)間步長(zhǎng)下耦合計(jì)算結(jié)果較為吻合，說明該耦合方法可用于管道瞬態(tài)過程模擬。

圖3 閥門測(cè)壓管水頭Figure 3 Valve gauge head

3 工程應(yīng)用

3.1 工程計(jì)算模型

以某抽水蓄能電站為例，該電站為一洞四機(jī)方案，設(shè)有下游調(diào)壓室，管道示意圖如圖4 所示。機(jī)組號(hào)為1 號(hào)、2 號(hào)、3 號(hào)和4 號(hào)，機(jī)組及管道其他參數(shù)如表1 所示。

表1 某抽水蓄能電站運(yùn)行參數(shù)Table 1 Operating parameters of a pumped storage power plant

表2 某抽水蓄能電站運(yùn)行參數(shù)Table 2 Operating parameters of a pumped storage power plant

圖4 某抽水蓄能電站管線布置示意圖Figure 4 Schematic diagram of pipeline layout of a pumped storage power plant

本文研究抽水蓄能電站尾水管空化和液柱分離受安裝高程的影響情況。尾水管部分采用FVM 方法建模計(jì)算，其他部分采用MOC 方法計(jì)算。同時(shí)，為了兼顧計(jì)算效率和準(zhǔn)確性，F(xiàn)VM 部分取0.001s 時(shí)間步長(zhǎng)計(jì)算，MOC 部分取0.01s 時(shí)間步長(zhǎng)計(jì)算。

計(jì)算工況為該工程過渡過程大波動(dòng)尾水管進(jìn)口壓強(qiáng)的最不利工況CT4-5：上庫正常蓄水位1262.5m，下庫死水位873m，四臺(tái)發(fā)電機(jī)組均以額定出力運(yùn)行。在0時(shí)刻，同一水力單元3號(hào)、4 號(hào)機(jī)組甩符荷，1 號(hào)、2 號(hào)機(jī)組5s 后甩負(fù)荷，導(dǎo)葉緊急關(guān)閉。

3.3 安裝高程敏感性分析

為研究空化和液柱分離對(duì)于抽水蓄能電站尾水管進(jìn)口壓強(qiáng)的影響，本次數(shù)值模擬，選取了6 組不同機(jī)組安裝高程進(jìn)行計(jì)算，并對(duì)1 號(hào)機(jī)組尾水管進(jìn)口壓強(qiáng)進(jìn)行分析。計(jì)算結(jié)果以相對(duì)壓強(qiáng)對(duì)應(yīng)的壓強(qiáng)水頭為準(zhǔn)。

6 組數(shù)值仿真模擬試驗(yàn)計(jì)算結(jié)果如圖5 所示。顯然，在不考慮壓力脈動(dòng)和計(jì)算誤差的情況下，安裝高程從初始的805m 上升到825m 時(shí)，尾水管進(jìn)口壓強(qiáng)隨著安裝高程上升而不斷下降，波形完全一致。當(dāng)機(jī)組安裝高程上升到830m 時(shí)，尾水管進(jìn)口壓強(qiáng)會(huì)下降到汽化壓強(qiáng)以下，進(jìn)而發(fā)生空化。

圖5 不同安裝高程尾水管進(jìn)口壓強(qiáng)水頭Figure 5 Tailpipe inlet pressure head at different installation elevations

3.3 水柱分離模擬對(duì)比

圖6 給出了初始安裝高程為830m 條件下，運(yùn)用FVM 采用空化模型和不采用空化模型兩種方法對(duì)1 號(hào)機(jī)組尾水管進(jìn)行數(shù)值模擬的計(jì)算結(jié)果。采用空化模型計(jì)算，尾水管進(jìn)口壓強(qiáng)達(dá)到汽化壓強(qiáng)時(shí)，壓強(qiáng)不再下降，出現(xiàn)空化和液柱分離，之后的最大升壓相比初始?jí)簭?qiáng)上升62.85m。而不發(fā)生空化的數(shù)學(xué)模型計(jì)算結(jié)果最大升壓相比初始?jí)簭?qiáng)上升38.15m，遠(yuǎn)小于發(fā)生液柱分離后的壓強(qiáng)升高值。并且發(fā)生液柱分離后壓強(qiáng)上升速率更快。同時(shí)，發(fā)生液柱分離后會(huì)使得波速降低，進(jìn)而水錘波動(dòng)周期向后推移。在導(dǎo)葉完全關(guān)閉后，尾水管進(jìn)口壓強(qiáng)趨于穩(wěn)定，與無空化模型計(jì)算結(jié)果相似。

圖6 尾水管進(jìn)口壓強(qiáng)水頭Figure 6 Tail pipe inlet pressure head

圖7 單獨(dú)給出了通過FVM 空化模型計(jì)算得到的尾水管進(jìn)口壓強(qiáng)隨時(shí)間變化曲線圖。在0 時(shí)刻，3 號(hào)、4 號(hào)機(jī)組導(dǎo)葉瞬時(shí)關(guān)閉，受此影響，機(jī)組尾水管壓強(qiáng)緩慢下降。5s 時(shí)刻，該機(jī)組導(dǎo)葉開始關(guān)閉，尾水管壓強(qiáng)開始迅速下降，并在9.35s 時(shí)刻達(dá)到汽化壓強(qiáng)，壓強(qiáng)停止下降，尾水管進(jìn)口發(fā)生空化，進(jìn)而出現(xiàn)液柱分離。10.93s 時(shí)刻，如圖7 中1 點(diǎn)所示，下游增壓波到達(dá)尾水管進(jìn)口，液柱彌合，同時(shí)液柱分離處出現(xiàn)兩個(gè)反射波，增壓波到達(dá)液柱分離處，受空穴反射為減壓波向下游傳播，而在這同時(shí)，增壓波繼續(xù)向機(jī)組傳播，在機(jī)組處反射為增壓波向下游傳播。受到增壓波影響，尾水管壓強(qiáng)迅速上升，如圖7 中2 點(diǎn)所示。在13.48s 時(shí)刻，二次反射增壓波到達(dá)尾水管進(jìn)口處，形成液柱分離水錘的二次升壓，升壓極值為95.21m，遠(yuǎn)大于不發(fā)生空化的升壓極值70.52m。

圖7 尾水管進(jìn)口壓強(qiáng)水頭Figure 7 Tail pipe inlet pressure head

4 結(jié)論

本文基于特征線法和有限體積法，引入相變空化模型并耦合邊界條件方程，提出了一種抽水蓄能機(jī)組水柱分離計(jì)算新方法，在不考慮壓力脈動(dòng)和計(jì)算誤差情況下，實(shí)現(xiàn)了對(duì)該電站尾水管空化過程的一維數(shù)值模擬，并闡明了水柱分離過程中尾水管進(jìn)口壓強(qiáng)的動(dòng)態(tài)演變規(guī)律。該方法克服了傳統(tǒng)特征線法不能模擬抽水蓄能水柱分離現(xiàn)象的局限性，為抽水蓄能電站安裝高程設(shè)計(jì)及安全穩(wěn)定運(yùn)行分析提供了有效工具。