基于一維N—S方程的池蘆河道型水庫(kù)調(diào)洪計(jì)算

□李雙江 □趙正鵬 □周 斌

(1三門(mén)峽市水利勘測(cè)設(shè)計(jì)有限責(zé)任公司；2四川省廣元水文水資源勘測(cè)局；3廣東省汕尾市水利水電規(guī)劃設(shè)計(jì)院)

1 基于一維N-S方程的調(diào)洪計(jì)算方法

近年來(lái)隨著計(jì)算能力的逐步提高，一維N—S方程在工程中的應(yīng)用也越來(lái)越普遍，積累了大量的應(yīng)用經(jīng)驗(yàn)，已成為一個(gè)可靠的工程分析工具。一維N—S基本方程為：

式中：z—水位；Q—流量；t—時(shí)間；b—河寬；x—斷面間距；g—重力加速度；A—過(guò)水面積；Sf—水力坡降。

一維N—S方程常采用差分法求解。文章采用Preissmann差分格式，該格式的表達(dá)式為：

其中系數(shù)和阻力項(xiàng)表達(dá)式為：

式中：j—計(jì)算節(jié)點(diǎn)號(hào)；n—計(jì)算時(shí)段；θ—變權(quán)系數(shù)；△x—分段長(zhǎng)度；△t—分段時(shí)長(zhǎng)；f—函數(shù)；

對(duì)于河道型水庫(kù)，可將河道型水庫(kù)的庫(kù)區(qū)劃為n個(gè)斷面，差分后一維N—S方程將形成有2n個(gè)變量、2n-2階非線性方程組，必須引入初始條件和兩個(gè)邊界條件，才可以求解。假定洪水發(fā)生前（t=0）各斷面為基流，水位為基流下的恒定流水面線，可作為求解的初始條件；上游（庫(kù)尾）邊界條件為已知時(shí)刻的入庫(kù)洪水流量，下游（壩址）邊界條件為已知水位—流量關(guān)系（泄流曲線）；由此可建立河道型水庫(kù)調(diào)洪的N—S非線性方程組，逐時(shí)段求解。

2 池蘆水庫(kù)采用一維N-S方程調(diào)洪計(jì)算過(guò)程

池蘆水庫(kù)位于黃河流域洛河水系支流永昌河上，壩址控制流域面積18.10 km2，水庫(kù)總庫(kù)容310.50萬(wàn)m3，是一座以灌溉、供水為主，兼有防洪任務(wù)的小（1）型水庫(kù)，防洪標(biāo)準(zhǔn)為30 a一遇設(shè)計(jì)、300 a一遇校核。池蘆水庫(kù)庫(kù)區(qū)兩岸水面最大寬度約150m，庫(kù)區(qū)長(zhǎng)約1.59 km，庫(kù)區(qū)內(nèi)無(wú)較大的支流匯入，是典型的河道型水庫(kù)。

池蘆水庫(kù)的一維N—S差分非線性方程組采用擬牛頓法迭代求解。將（1）式、（2）式采用Preissmann隱式差分法并擬牛頓化，可改寫(xiě)為：

式中：j—計(jì)算節(jié)點(diǎn)號(hào)；n—計(jì)算時(shí)段；K—流量模數(shù)。

將池蘆水庫(kù)長(zhǎng)1.59 km的庫(kù)區(qū)劃分為16個(gè)斷面（壩前樁號(hào)為0+000，庫(kù)尾樁號(hào)為1+592）；初始狀態(tài)取各斷面水位為正常蓄水位616m，流量為0m3/s；上游（庫(kù)尾）為已知流量（設(shè)計(jì)洪水過(guò)程線）邊界；下游（壩址）為已知水位—流量關(guān)系（池蘆水庫(kù)溢洪道的泄流曲線）邊界。取糙率n=0.034，△t=15 s，將16個(gè)斷面數(shù)據(jù)、初始條件和邊界條件代入（6）、（7）式采用追趕法迭代求解，經(jīng)逐時(shí)段演算，即可解得整個(gè)調(diào)洪過(guò)程。

3 調(diào)洪成果及對(duì)比

一維N—S方程求解的調(diào)洪成果和水量平衡法的調(diào)洪成果見(jiàn)表1，一維N—S方程求解成果的校核洪水工況的水位、流量、流速過(guò)程見(jiàn)圖1。兩種方法計(jì)算校核洪水工況的入庫(kù)洪水、壩前水位、泄量過(guò)程見(jiàn)圖2。

表1 池蘆水庫(kù)調(diào)洪成果表

采用N—S方程調(diào)洪的水位過(guò)程與水量平衡法的調(diào)洪水位成果極為接近，兩者的壩前水位過(guò)程線和泄流過(guò)程線基本重合；N—S方程計(jì)算的庫(kù)區(qū)水位基本水平，壩前（0+000）和庫(kù)尾（1+592）的水位基本同步，反映出水庫(kù)水深相對(duì)較大、庫(kù)區(qū)底坡較陡，池蘆水庫(kù)的動(dòng)庫(kù)容效應(yīng)不明顯；庫(kù)尾段流速變化大，校核工況庫(kù)尾的流速為1.72m/s，壩前僅為0.13m/s，相應(yīng)庫(kù)尾的流速水頭達(dá)0.15m，壩前卻<1.00×10-3m，致使洪峰段庫(kù)尾的水位出現(xiàn)了輕微的下凹現(xiàn)象。

4 一維N—S方程法和水量平衡方程法運(yùn)用的差別

圖1 校核工況調(diào)洪的水位、流量、流速變化過(guò)程圖

圖2 校核洪水工況入庫(kù)洪水、壩前水位、泄量變化過(guò)程圖

從池蘆水庫(kù)調(diào)洪的成果來(lái)看，水量平衡法和N—S方程法出入不大，精度總體上相等。通過(guò)池蘆水庫(kù)對(duì)兩種調(diào)洪方法的使用，揭示出兩種方法存在較大的不同，計(jì)算所需的資料、適應(yīng)的水庫(kù)類型、對(duì)動(dòng)庫(kù)容的使用和可獲得的計(jì)算成果上都是不相同的。兩種方法的主要區(qū)別在于：水量平衡法計(jì)算需要庫(kù)容曲線，適用于任意類型水庫(kù)，通常采用靜態(tài)庫(kù)容，計(jì)算成果為壩前水位過(guò)程和水庫(kù)洪水過(guò)程；而N—S方程法計(jì)算需要河道橫斷，適用于河道型水庫(kù)，可自動(dòng)計(jì)入動(dòng)庫(kù)容作用，能夠顯示出庫(kù)洪水過(guò)程及整個(gè)庫(kù)區(qū)內(nèi)的水位、流量及流速變化過(guò)程。

5 結(jié)語(yǔ)

目前工程設(shè)計(jì)中一般采用水量平衡法進(jìn)行水庫(kù)調(diào)洪計(jì)算，通?？梢詽M足工程需要；但一維N—S方程也是一種優(yōu)良的調(diào)洪算法，雖然目前N—S方程法僅能適用于河道型水庫(kù)調(diào)洪，但能自動(dòng)計(jì)入動(dòng)庫(kù)容作用，能揭示更為豐富的水力細(xì)節(jié)，也是一種可資使用的調(diào)洪方法。

[1]巨江.工程水力學(xué)數(shù)值仿真與可視化[M].北京：中國(guó)水利水電出版社，2010.

[2]汪德?tīng)?計(jì)算水力學(xué)理論與應(yīng)用[M].北京：科學(xué)出版社，2011.

[3]鄭咸義，姚仰新，雷秀仁，等.應(yīng)用數(shù)值分析[M].廣州：華南理工大學(xué)出版社，2008.

[4]周斌，張彥軍.納潮對(duì)龍津河水體更新的模擬分析[J].人民珠江,2015（3）:23-25.

[5]李雙江，楊莉，周建波，等.陜縣池蘆抗旱小水庫(kù)水源工程項(xiàng)目建議書(shū)[R].三門(mén)峽：三門(mén)峽市水利勘測(cè)設(shè)計(jì)有限責(zé)任公司，2015.