李甲振，王 濤，郭永鑫，薛興祖，紀(jì)昌知

（1. 中國水利水電科學(xué)研究院 流域水循環(huán)模擬與調(diào)控國家重點(diǎn)實(shí)驗(yàn)室，北京 100038； 2. 吉林省水利水電勘測設(shè)計(jì)研究院，吉林 長春 130021； 3. 中國三峽建工（集團(tuán)）有限公司，浙江 杭州 310052）

0 引 言

長距離調(diào)水是優(yōu)化水資源空間配置、緩解局部地區(qū)或農(nóng)業(yè)灌溉水資源短缺的工程措施，輸水方式一般有管、涵、隧有壓和渠、隧無壓兩種［1］。一般情況下，調(diào)水工程沿線會(huì)設(shè)置若干分水口，通過泵站加壓或重力自流方式向城鎮(zhèn)、工農(nóng)業(yè)供水，如引漢濟(jì)渭［2］、滇中引水［3，4］、大伙房水庫輸水［5］等工程。當(dāng)計(jì)算分水口的水力過渡過程時(shí)，多數(shù)工程需對(duì)干線和支線同時(shí)求解，這就會(huì)碰到有壓、無壓輸水耦合的問題。

針對(duì)有壓、無壓耦合輸水系統(tǒng)的水力過渡過程，一種常見的求解方法是采用明渠非恒定流方程和窄縫法描述管道的有壓流。楊開林［6］采用Preissmann 四點(diǎn)隱式差分算法進(jìn)行求解，并給出了保證計(jì)算收斂的措施和方法；劉梅清等［7］則采用特征線方法進(jìn)行求解，萬五一等［8］給出了最優(yōu)時(shí)間步長和管、渠空間步長的確定方法。王衍超等［9］在處理南水北調(diào)北京段的明滿流耦合問題時(shí)，也是采用了窄縫法進(jìn)行計(jì)算。莫鐵祥和李國棟［10］針對(duì)管流負(fù)壓工況，提出了一種結(jié)合CFL 約束條件的流態(tài)判別方式。另一種求解方法則是分別計(jì)算有壓、無壓耦合系統(tǒng)的水力過渡過程，連接處的邊界條件進(jìn)行特殊處理。李占松等［11］區(qū)分管渠、渠管兩種情形，給出了斷面結(jié)合處特征線法內(nèi)邊界的處理方式。

對(duì)于有壓、無壓耦合的長距離調(diào)水工程，無壓輸水系統(tǒng)的計(jì)算時(shí)長一般為數(shù)小時(shí)至數(shù)天，時(shí)間步長一般為數(shù)秒至數(shù)分；有壓輸水系統(tǒng)考慮特征線方法的計(jì)算穩(wěn)定性和計(jì)算精度，時(shí)間步長一般小于0.1 s，甚至是0.001 s。如果兩者采用同一時(shí)間步長，會(huì)占用大量的CPU，降低計(jì)算效率。因此，有必要對(duì)無壓輸水系統(tǒng)和有壓輸水系統(tǒng)采用不同的時(shí)間步長進(jìn)行求解。

本文以某長距離調(diào)水工程為例，研究給出無壓、有壓輸水系統(tǒng)耦合的求解算法，分析泵站事故斷電、正常啟動(dòng)、停機(jī)的控制策略和水力特性，供類似工程的水力過渡過程計(jì)算參考。

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

某調(diào)水工程從上游水庫引水，通過長224 km的無壓隧洞輸送至下游水庫。無壓隧洞底坡為1/4010，斷面為馬蹄形，尺寸為6.56 m×6.56 m，設(shè)計(jì)流量為50.0 m3/s。在輸水隧洞126 km處，向下開挖穩(wěn)壓調(diào)節(jié)池，側(cè)向引水進(jìn)入泵站，設(shè)計(jì)流量為10.2 m3/s。隧洞底高程為242.20 m，穩(wěn)壓調(diào)節(jié)池深6.44 m，兩側(cè)邊坡1∶10，底長30.00 m，頂長158.80 m，與明流隧洞同寬6.56 m［圖1（a）］。

圖1 工程布置圖（單位：m）Fig.1 Structure layout of pump station diversion

側(cè)向引水通過有壓隧洞進(jìn)入閘門井，隧洞長273.00 m，底坡為1/54；斷面為圓形，直徑為3.00 m。閘門井的順?biāo)鏖L度為8.00 m，寬度為3.00 m。閘門井后接集管進(jìn)行分流，兩側(cè)對(duì)稱布置4 根進(jìn)水管，長度分別為45.00 m、31.00 m；斷面為圓形，直徑為2.00 m［圖1（b）］。

該泵站引水系統(tǒng)包括了無壓輸水和有壓輸水，控制方程和求解算法如下。

1.1 控制方程

無壓輸水系統(tǒng)的控制方程為圣維南方程組，包括連續(xù)性方程和動(dòng)量方程：

式中：A為過流面積，m2；t為時(shí)間變量，s；Q為流量，m3/s；x為空間變量，m；q為單位渠道長度的側(cè)向流量，m3/s；β為斷面流速分布不均引入的修正系數(shù)；g為重力加速度，m2/s；h為水深，m；S0為河床底坡；Sf為摩阻比降，計(jì)算公式為：

式中：K為流量模數(shù)，；n為糙率；R為水力半徑，m。

有壓輸水系統(tǒng)的控制方程包括運(yùn)動(dòng)方程和連續(xù)方程：

式中：H為從基準(zhǔn)線算起的測壓管水頭，m；V為管道中平均流速，m/s；f為Darcy-Weisbach 沿程水力損失系數(shù)；D為管道直徑，m；a為水錘波速，m/s；θ為管軸傾角。

1.2 求解算法

無壓輸水系統(tǒng)的控制方程組，常采用Preissmann 四點(diǎn)隱式差分法格式求解；有壓輸水系統(tǒng)的控制方程組，常采用特征線方法進(jìn)行求解。文獻(xiàn)［12-14］均對(duì)離散方法、計(jì)算程序以及基本的水庫、調(diào)壓塔等邊界條件進(jìn)行了詳細(xì)介紹，此處不再贅述。數(shù)值計(jì)算的難點(diǎn)在于，如何聯(lián)立求解無壓輸水系統(tǒng)和有壓輸水系統(tǒng)的水力過渡過程。

在無壓輸水系統(tǒng)的一個(gè)時(shí)間步長ΔT內(nèi)，穩(wěn)壓調(diào)節(jié)池水深的改變量Δh相對(duì)于水泵揚(yáng)程是一個(gè)小量，對(duì)有壓輸水系統(tǒng)的水力過渡過程影響很小，可忽略不計(jì)。但管道系統(tǒng)與穩(wěn)壓調(diào)節(jié)池之間的水量交換，對(duì)明流隧洞的水力過渡過程影響較大，必須考慮。因此，本文采用了一種弱耦合的方式進(jìn)行求解，即，①計(jì)算有壓輸水系統(tǒng)T～T+ ΔT的水力過渡過程，計(jì)算步數(shù)為N=ΔT/Δt，穩(wěn)壓調(diào)節(jié)池水位恒定為時(shí)間T的數(shù)值。② 計(jì)算N個(gè)Δt時(shí)間步長，有壓隧洞首斷面的流入/流出的流量。將N個(gè)流量的平均值作為該節(jié)點(diǎn)（岔點(diǎn)）的邊界條件，求解T+ ΔT時(shí)刻的水位和流量。③ 轉(zhuǎn)入步驟①、②繼續(xù)求解。

數(shù)值求解的計(jì)算流程如圖2所示。

圖2 計(jì)算流程圖Fig.2 Standard calculation procedure

1.3 典型建筑物處理問題

泵站從無壓輸水隧洞的中間部分引水，通過有壓隧洞、閘門井、集管、進(jìn)水管后進(jìn)入水泵。工程中的閘門前后均為有壓隧洞，規(guī)劃的功能是截?cái)嗨?，正常運(yùn)行時(shí)處于全開位置。當(dāng)發(fā)生事故斷電、水泵正常啟動(dòng)或停運(yùn)時(shí)，由于閘門井具有自由水面，起到了調(diào)壓井的作用。因此，閘門井需作為雙向調(diào)壓塔考慮。當(dāng)水位高于244.60 m 時(shí)，調(diào)壓塔水面面積為3.00 m×8.00 m；當(dāng)水位為239.60～244.60 m時(shí)，調(diào)壓塔水面面積線性變化。

1.4 模型校驗(yàn)

青草沙水庫的原水通過2 根長13 650 m、直徑5.5 m 的輸水隧洞輸送至五號(hào)溝泵站前池，繼而向下游各用戶配水。某次水泵調(diào)試運(yùn)行時(shí)，李靜毅等［15］觀察到配水池異常的水位波動(dòng)，振幅高達(dá)1m。本文對(duì)該過程進(jìn)行了仿真計(jì)算，關(guān)鍵參數(shù)如表1所示。仿真計(jì)算的水位波動(dòng)、峰谷值和周期與現(xiàn)場實(shí)測結(jié)果一致，如圖3所示，說明模型能夠較好地模擬輸水系統(tǒng)的水力過渡過程特性。

表1 算例參數(shù)Tab.1 Parameters of case study

圖3 仿真計(jì)算與現(xiàn)場實(shí)測對(duì)比Fig.3 Comparisons between numerical simulation and field test

2 結(jié)果與討論

2.1 水泵事故斷電的水力過渡過程

設(shè)計(jì)工況下，穩(wěn)壓調(diào)節(jié)池的上游流量為50.0 m3/s，下游流量為39.8 m3/s，泵站的取水流量為10.2 m3/s。發(fā)生事故斷電后，水泵出口閥門兩段關(guān)閥，第一段快關(guān)15 s，第二段時(shí)間45 s，拐點(diǎn)20%全開度；上游水庫的放水流量在100 s內(nèi)由50.0 m3/s減小為39.8 m3/s。

閘門井水位波動(dòng)如圖4 所示。泵站正常運(yùn)行時(shí)，閘門井水位為246.20 m。事故斷電后，水泵出口的閥門兩階段關(guān)閉，閘門井水位迅速增加，32 s 達(dá)到最高值247.65 m，水位升高0.45 m；之后閘門井水位降低，82 s 達(dá)到極小值246.42 m。0～40 h 內(nèi)，閘門井水位在0～9 h 內(nèi)振蕩上升，9～30 h 內(nèi)逐漸降低，之后穩(wěn)定為246.40 m。主要原因是，水庫下泄50.0 m3/s、泵站取用10.2 m3/s時(shí)，上游126 km 隧洞基本為50.0 m3/s 對(duì)應(yīng)的設(shè)計(jì)水深4.75 m，下游98 km 隧洞基本為39.8 m3/s 對(duì)應(yīng)的設(shè)計(jì)水深3.99 m。水泵事故斷電、閥門關(guān)閉后，泵站停止取水，上游126 km隧洞內(nèi)的水向穩(wěn)壓調(diào)節(jié)池及下游輸送，致使穩(wěn)壓調(diào)節(jié)池、閘門井的水位逐漸增加；9 h 后，上游水庫下泄水量減小產(chǎn)生的減壓波傳播至穩(wěn)壓調(diào)節(jié)池，閘門井水位開始降低。

圖4 閘門井水位波動(dòng)Fig.4 Water level variation of gate shaft

穩(wěn)壓調(diào)節(jié)池上、下游明流隧洞的流量過程如圖5 所示。水泵事故斷電后，上游水庫的放水流量在100 s 內(nèi)減小，隧洞20、40、60、80 和100 km 處的流量逐漸減小，水位逐漸降低，約22 h后，100 km 處的流量達(dá)到40.2 m3/s［圖5（a）］。由于水泵出口閥門關(guān)閉，側(cè)向引水減小為零，穩(wěn)壓調(diào)節(jié)池上游隧洞3.99～4.75 m（39.8、50.0 m3/s 對(duì)應(yīng)的設(shè)計(jì)水深）之間的水體只能輸送至下游，致使下游隧洞的流量逐漸增加。當(dāng)上游水庫下泄水量減小產(chǎn)生的減壓波傳播至相應(yīng)位置時(shí)，流量開始減?。蹐D5（b）］。穩(wěn)壓調(diào)節(jié)池、下游20 km、下游40 km、下游60 km 和出口的流量最大值分別為49.6、49.3、48.9、48.6 和48.2 m3/s。相對(duì)比設(shè)計(jì)流量39.8 m3/s，40 h內(nèi)多輸送至下游水庫的流量為55.2 萬m3。

圖5 明流隧洞流量波動(dòng)Fig.5 Discharge variation of free-flow tunnel

穩(wěn)壓調(diào)節(jié)池上、下游明流隧洞的水位波動(dòng)如圖6 所示。水泵事故斷電后，隨著閥門兩階段關(guān)閉，穩(wěn)壓調(diào)節(jié)池的水位在36 s（0.01 h）迅速抬升，72 s（0.02 h）升高0.17 m［圖6（a）］。穩(wěn)壓調(diào)節(jié)池側(cè)向分水減小引起的水位增高，致使上游200 m 處水位略降低后增加。與閘門井水位波動(dòng)類似，0～40 h 內(nèi)，穩(wěn)壓調(diào)節(jié)池水位在0～9 h 內(nèi)振蕩上升，9～30 h 內(nèi)逐漸降低，之后穩(wěn)定［圖6（b）］。數(shù)值計(jì)算過程中，一個(gè)時(shí)間步長ΔT=2 s內(nèi)穩(wěn)壓調(diào)節(jié)池的最大水位變幅為0.02 m，對(duì)有壓隧洞水力過渡過程的影響可忽略，這也說明了本文提出的數(shù)值計(jì)算方法是可行的。

2.2 水泵正常啟停的水力過渡過程

水泵正常啟動(dòng)時(shí)，輸水系統(tǒng)的初始狀態(tài)為：穩(wěn)壓調(diào)節(jié)池的上游流量、下游流量均為39.8 m3/s；輸水系統(tǒng)的最終狀態(tài)為：穩(wěn)壓調(diào)節(jié)池的上游流量為50.0 m3/s，下游流量為39.8 m3/s，泵站的取水流量為10.2 m3/s。經(jīng)研究給出的控制策略為：上游水庫的放水流量在100 s內(nèi)由39.8 m3/s 調(diào)整為50.0 m3/s；三臺(tái)機(jī)組的泵后閥門在20 h開始動(dòng)作，100 s內(nèi)線性開啟。

閥門未開啟時(shí)，閘門井水位為246.95 m。閥門開啟后，水泵吸水產(chǎn)生的壓力波動(dòng)在穩(wěn)壓調(diào)節(jié)池、閘門井之間傳播，使閘門井水位產(chǎn)生一定的波動(dòng)，最大值為247.11 m，最小值為246.00 m，振幅為1.11 m［圖7（a）］。泵后閥門開啟所產(chǎn)生的閘門井水位波動(dòng)，15 min 后振幅趨于0。0～20 h 內(nèi)，上游放水增加使閘門井水位增加；20～25 h，閘門井水位逐漸降低，之后趨于穩(wěn)定水位246.20 m［圖7（b）］。

圖7 閘門井水位波動(dòng)Fig.7 Water level variation of gate shaft

上游水庫的放水流量增加，隧洞20、40、60、80 和100 km 處的流量逐漸增大，水位逐漸升高，約22 h后，100 km 處的流量達(dá)到49.5 m3/s［圖8（a）］。閥門開啟后，水泵取水，穩(wěn)壓調(diào)節(jié)池流量迅速減小至設(shè)計(jì)流量39.8 m3/s。下游20 km、下游40 km、下游60 km 和出口的流量最大值分別為46.9、45.8、45.0 和44.3 m3/s。37 h 后，下游出口流量恢復(fù)為40.2 m3/s［圖8（b）］。相對(duì)比設(shè)計(jì)流量39.8 m3/s，40 h內(nèi)多輸送至下游水庫的流量為17.61 萬m3。

圖8 明流隧洞流量波動(dòng)Fig.8 Discharge variation of free-flow tunnel

水泵正常停機(jī)時(shí)，經(jīng)研究給出的控制策略為：上游水庫的放水流量在100 s內(nèi)由50.0 m3/s 調(diào)整為39.8 m3/s；三臺(tái)機(jī)組的泵后閥門在20 h 后開始動(dòng)作，第一段快關(guān)15 s，第二段時(shí)間45 s，拐點(diǎn)20%全開度。其水力過渡過程與正常啟動(dòng)類似，不再贅述。

3 結(jié) 論

針對(duì)某泵站從長距離明流隧洞中間部分引水的水力過渡過程計(jì)算問題，提出了一種弱耦合的求解方式，即，認(rèn)為連接處的水位在ΔT（無壓輸水系統(tǒng)的時(shí)間步長）內(nèi)是恒定的，計(jì)算有壓輸水系統(tǒng)N個(gè)時(shí)間步（N=ΔT/Δt，Δt為有壓輸水系統(tǒng)的時(shí)間步長）的水力瞬變過程，有壓輸水系統(tǒng)的進(jìn)出流量影響無壓輸水系統(tǒng)。給出了水泵事故斷電、正常啟動(dòng)和正常停機(jī)工況下，泵后閥門、上游水庫的控制策略及其對(duì)下游水庫的影響。

上述工程中，閘門井具有自由水面，前后為長度較短的有壓輸水隧洞（273.00 m）和進(jìn)水管（45.00 m、31.00 m），在水力過渡過程中起到了調(diào)壓塔的功能。雖然命名閘門井，但不能僅作為截流建筑，需具體問題具體分析，還應(yīng)考慮其調(diào)壓作用。