江穎 李珍 安靜泊

摘?要：小浪底水利樞紐進水塔群前的泥沙淤積影響工程的安全運行和水庫的發(fā)電效益，制定合理的泄水孔洞調(diào)度方案至關(guān)重要。為此采用立面二維水沙數(shù)學(xué)模型開展小浪底水庫底孔防淤堵的數(shù)值模擬研究。在設(shè)計水沙條件下，計算和比較兩種泄水孔洞運用方案對壩區(qū)河床縱剖面淤積形態(tài)、排沙比和發(fā)電效益的影響，結(jié)果表明：兩種方案均能滿足進水塔群防淤堵的要求，為盡可能發(fā)揮小浪底水庫的發(fā)電效益，泄水孔洞調(diào)度方案推薦方案1，即出庫流量小于發(fā)電洞泄流量時優(yōu)先啟用發(fā)電洞泄流，當(dāng)出庫流量超過發(fā)電洞泄流量時超出部分盡量通過排沙洞、明流洞、孔板洞泄流，減少庫區(qū)淤積。

關(guān)鍵詞：小浪底水庫;防泥沙淤堵;立面二維水沙數(shù)學(xué)模型;泄水孔洞調(diào)度方案

中圖分類號：TV697.1;TV882.1?文獻標(biāo)志碼：A

doi：10.3969/j.issn.1000-1379.2020.10.003

Comparison of Operation Dispatch Scheme of the Discharge Outlets in the Xiaolangdi Reservoir

JIANG Ying， LI Zhen， AN Jingbo

（Xiaolangdi Multi-Purpose Dam Project Management Center， Ministry of Water Resources， Zhengzhou 450000， China）

Abstract：Preventing sediment siltation in front of intake tower group of the Xiaolangdi Reservoir affects the safety operation and the generation benefits of the project. It is vital to determine a proper operation dispatch scheme of the outlet structures. This paper did the numerical research on the prevention of silt plug in the Xiaolangdi Reservoir by a vertical 2D numerical model. The bed evolution， sediment delivery ratio and generation benefit were discussed under typical water and sediment conditions with the two different dispatch schemes of the outlets. The results of the numerical modeling show that neither of the dispatch schemes will cause the sediment siltation and the scheme 1 shows a better performance in promoting the generation benefits. The present work provides the technical support for the operation dispatch scheme of outlets in engineering practices.

Key words： Xiaolangdi Reservoir; preventing sediment siltation; vertical 2D numerical model; operation dispatch scheme of outlet structures

黃河小浪底水利樞紐是黃河干流的關(guān)鍵控制性骨干工程，在黃河治理開發(fā)中具有十分重要的戰(zhàn)略地位。為了滿足以防洪（包括防凌）、減淤為主，兼顧供水、灌溉、發(fā)電等的開發(fā)任務(wù)需要，小浪底水利樞紐共設(shè)置了3條孔板洞、3條排沙洞、3條明流洞、6條發(fā)電洞和1條灌溉洞[1]，各隧洞設(shè)有工作、事故、檢修閘門，控制各條泄水洞的進水“咽喉”，形成底部泄洪排沙、中間引水發(fā)電和上部泄洪排漂的格局[2]，擔(dān)負著完成樞紐任務(wù)和確保工程安全的重要使命。

黃河來水含沙量高且來沙量大。小浪底水庫運行至2019年10月，庫區(qū)三角洲頂點高程已經(jīng)達到212.70 m，頂點距壩里程約10.32 km。壩前1.32 km處淤積高程由137.50 m抬高至181.64 m，淤積抬升了44.14 m，高于最低進水口底板高程6.64 m;進水塔前60 m處淤積高程為176.18 m，高于最低進水口底板高程1.18 m。受泥沙淤積影響，有可能導(dǎo)致進水口淤堵、閘門啟閉困難等。如2019年7月17日，在執(zhí)行全關(guān)2號排沙洞事故閘門操作時，發(fā)現(xiàn)閘門底坎有厚約330 mm的淤積泥沙，致使無法全關(guān)到位。據(jù)統(tǒng)計，2018—2019年小浪底工程事故閘門啟閉遇到的淤堵及卡阻問題達12次之多。

隨著小浪底水庫庫區(qū)泥沙淤積發(fā)展，三角洲頂點進一步向壩前推進，淤積面不斷抬高，泄水建筑物進水口泥沙淤堵的風(fēng)險也越來越大，影響工程效益發(fā)揮和樞紐工程安全運行，開展小浪底水利樞紐進水塔群前防淤堵研究十分必要。目前已開展了一系列的模型試驗研究[3-6]，并取得了一定的研究成果，但是物理模型試驗的周期長，且對不同工況的靈活度較低。筆者利用立面二維數(shù)學(xué)模型，在典型設(shè)計水沙條件下對不同泄水孔洞運用方案的底孔防淤堵進行數(shù)值模擬研究，從河床縱剖面、排沙比、發(fā)電效益等方面綜合分析比較，從而選定合理的泄水孔洞運用方案，為確保水庫安全運行并充分發(fā)揮小浪底水庫的工程效益提供技術(shù)支持。

1?立面二維數(shù)學(xué)模型及其控制方程

近年來三維水流數(shù)學(xué)模型發(fā)展迅速，可以直觀準確地描述水流的細節(jié)運動[7-9]，但是其計算成本非常高，很難普遍應(yīng)用于大尺度長歷時的實際工程問題研究。一維水流數(shù)學(xué)模型適用于長歷時長河段的計算，但其更側(cè)重于物理量在各斷面上的平均值，在某些特定的問題中難以反映物理要素在斷面上的分布。綜合而言，二維數(shù)學(xué)模型能夠較好地平衡計算量和水流細節(jié)描述的問題。二維模型分為平面二維和立面二維模型，對于水庫壩前段，伴隨著水庫泄水等調(diào)節(jié)控制運用，常有異重流發(fā)生，并且在孔口前易形成沖刷漏斗，不同深度的含沙量及流速差異明顯。為進行小浪底水利樞紐進水塔群前的防淤堵研究和泄水孔洞的調(diào)度方案比選，宜采用立面二維數(shù)學(xué)模型進行模擬和計算研究。

根據(jù)流體力學(xué)的基本守恒定律，分別以水流、泥沙以及河床為研究對象建立各自的控制方程[10]。

水流連續(xù)方程：

ux+vy=0（1）

水流運動方程：

ut+uux+vvy=-1ρpx+νt（2ux2+2uy2）（2）

vt+uvx+vvy=-g-1ρpy+νt（2vx2+2vy2）（3）

懸移質(zhì)運動方程：

hst+husx+h（v-ω）sy=ε（2hsx2+2hsy2）（4）

河床變形方程：

Zbt=Pγω（sa-s*a）γ′（5）

式中：x為沿水流方向坐標(biāo);y為垂直方向坐標(biāo)，向上為正方向;u為沿水流方向流速;v為垂直方向流速;p為壓強;ρ為清水密度;νt為紊流黏滯系數(shù);h為水深;ε為泥沙紊動擴散系數(shù);Zb為河床高程;ω為泥沙沉速;Pγ為泥沙孔隙率;s為含沙量;sa、s*a為床面含沙量及挾沙力;γ′為泥沙干容重。

為了使表示自由水面的水位函數(shù)ξ（x，t）和動量方程緊密聯(lián)系起來，將壓強分解為動水壓強和靜水壓強之和，即p=pd+ps。pd為動水壓強，是水體流動時因流線彎曲和流速不均勻所產(chǎn)生的附加壓強;ps為靜水壓強，ps=ρg[ξ（x，t）-y]，將壓強分解后，壓強梯度可分別改為

-px=-pdx-ρgξ（x，t）x（6）

-ρg-py=ρdy （7）人 民 黃 河?2020年第10期

2?數(shù)學(xué)模型的求解

2.1?方程離散

水流運動方程離散采用基于交錯網(wǎng)格離散的SIMPLE算法，如圖1所示，標(biāo)量變量和各種常數(shù)儲存在網(wǎng)格節(jié)點上，矢量變量及速度變量儲存在控制邊界上。

水流運動方程可離散為

un+1i+1，j+1/2=uni+1，j+1/2-（uux+vuy）Δt+

[-1ρpx+νt（2ux2+2uy2）]Δt（8）

vn+1i+1/2，j+1=vni+1/2，j+1-（uvx+vvy）Δt+

[-1ρpy+νt（2vx2+2vy2）]Δt（9）

式中：Δt為時間步長;n為時間節(jié)點編號;i為x方向節(jié)點編號;j為y方向節(jié)點編號。

懸移質(zhì)運動方程在矩形網(wǎng)格上采用有限體積法離散：

（hsk）n+1i+1/2，j+1/2-（hsk）ni+1/2，j+1/2Δt+

（hu）n+1i+1，j+1/2snki+1，j+1/2-（hu）ni+1/2，j+1/2snki，j+1/2Δx+（hv-hω）n+1i+1，j+1/2snki+1/2，j+1/2-（hu-hω）ni+1/2，j+1/2snki，j+1/2Δy=

CFFx-CFBxΔx+CFFy-CFByΔy（10）

式中：CFFx、CFBx、CFFy、CFBy為式（4）右端項的離散代數(shù)關(guān)系式，根據(jù)文獻[11]進行求解。

CFFx=εhn+1i+3/2，j+1/2+hn+1i+1/2，j+1/22sni+3/2，j+1/2-sni+1/2，j+1/2ΔxCFBx=εsni+3/2，j+1/2+sni+1/2，j+1/22hn+1i+3/2，j+1/2-hn+1i+1/2，j+1/2ΔxCFFy=εhn+1i+1/2，j+3/2+hn+1i+1/2，j+1/22sni+1/2，j+3/2-sni+1/2，j+1/2ΔxCFBy=εsni+1/2，j+3/2+sni+1/2，j+1/22hn+1i+1/2，j+3/2-hn+1i+1/2，j+1/2Δx（11）

河床變形方程離散形式為

ΔZb=Pγω（sa-s*a）γ′Δt（12）

式中：ΔZb為河床沖淤厚度。

2.2?封閉模式

泥沙運動極其復(fù)雜，對于其中河底挾沙力公式需要采用經(jīng)驗公式來進行封閉。立面二維水沙數(shù)學(xué)模型計算中與河床沖淤變形相關(guān)的因素是河床底部含沙量與挾沙力，通常選一接近河底的點作為參考點，參考點距離河底為ξa=（0.01～0.05）h，河底挾沙力公式為

sb=s*（1+gck）j1-gckj2（1ξ-1）z1（13）

式中：s*為挾沙力;c為謝才系數(shù)，c=h1/6/n;j1、j2分別為計算參數(shù)，j1=∫1ξa[（1-ξ）/ξ]Zdξ，j2=∫1ξa[（1-ξ）/ξ]Zln ξdξ，其中Z為懸浮指標(biāo)、ξ為水深;k為卡門常數(shù);z1=ωβkU*=Zβ，其中β為參數(shù)、U*為摩阻流速。

2.3?模型計算網(wǎng)格

模型計算網(wǎng)格布置如圖2所示。網(wǎng)格長度為50 m，高度為2 m，總數(shù)為90×50。

3?不同泄水孔洞調(diào)度方案比選

數(shù)值模擬的研究范圍為小浪底水庫庫區(qū)大峪河口至大壩段，在距離進水塔4.5 km范圍內(nèi)設(shè)置了31個漏斗測驗斷面，斷面編號自下而上依次為HHLD01—HHLD31，距塔210 m以內(nèi)斷面間距為50 m，距塔210～710 m范圍內(nèi)斷面間距為100 m，距塔710 m以上斷面間距為100～200 m。

3.1?典型設(shè)計水沙過程選取

20世紀80年代中后期以來，受人類活動和氣候變化的影響，黃河來水來沙條件發(fā)生了較大變化，尤其是2000年以來水沙條件變化更為明顯，來水來沙量大幅度減少。本次研究在1987年以來實測入庫水沙過程的基礎(chǔ)上選取入庫水沙條件，從不利于進水塔前防淤堵的角度考慮，選取平水偏豐沙過程。

由于小浪底水庫沖淤變化主要發(fā)生在汛期，因此針對7—9月對進水塔群前進行防淤堵模擬研究。通過對1987年以來7—9月三門峽站（入庫控制站）的水沙過程分析，選取平水偏豐沙的2003年7—9月實測入庫水沙過程作為入庫的典型水沙條件。

為了提高模擬精度，立面二維模型常采用較為精細的計算網(wǎng)格，但是受計算能力限制，模擬范圍往往局限于壩前區(qū)域。鑒于此，需要利用庫區(qū)一維水沙數(shù)學(xué)模型提供的壩前4.5 km處的流量、輸沙率過程，以此確定壩區(qū)立面二維水沙模型進口邊界的流量和含沙量分布。

由于沒有2003年的初始地形數(shù)據(jù)，因此基于2014年4月的初始地形資料，把2003年7—9月實測水沙過程作為計算條件，將計算結(jié)束后的地形作為庫區(qū)立面二維水沙數(shù)學(xué)模型計算的河床初始邊界條件，見表1和圖3。其次，把2003年7—9月進入壩區(qū)的實測水沙過程作為立面二維模型所采用的設(shè)計水沙過程，其設(shè)計流量過程如圖4所示，設(shè)計含沙量過程（該設(shè)計水沙過程以下簡稱為2003典型）見圖5，設(shè)計水沙量見表2。進水塔前高程與底孔進口高程均為175 m。

3.2?泄水孔洞調(diào)度方案

底孔防淤堵研究擬對兩種泄水孔洞調(diào)度方案進行比選。方案1：當(dāng)出庫流量小于發(fā)電洞泄流量時，優(yōu)先啟用發(fā)電洞泄流;當(dāng)出庫流量超過發(fā)電洞泄流量時，超出部分盡量通過排沙洞、明流洞、孔板洞泄流，減少庫區(qū)淤積。方案2：出庫流量大于400 m3/s時，發(fā)電洞過水80%，排沙洞分流20%。具體調(diào)度方案見表3。

3.3?數(shù)值模擬結(jié)果分析

采用庫區(qū)立面二維水沙數(shù)學(xué)模型，分別對2003典型水沙條件下泄水孔洞調(diào)度方案1和方案2進行計算和分析比較，提出推薦的泄水孔洞調(diào)度方案。

3.3.1?河床縱剖面分析

兩種泄水孔洞調(diào)度方案的河床淤積厚度和淤積量計算結(jié)果見表4，2003典型方案1和方案2河床縱剖面沖淤變化對比如圖6所示?？傮w而言，兩個方案計算結(jié)果差別不大，方案1和方案2在塔前約600 m范圍均形成沖刷漏斗，塔前60 m處與初始時刻（7月1日）相比淤積厚度分別為1.40 m和0.42 m，塔前4 500 m處與初始時刻相比淤積厚度分別為11.41 m和11.29 m。7月1日—9月30日，整個計算范圍平均淤積厚度分別為8.87 m和7.96 m，方案1與方案2相比，平均淤積厚度多0.91 m，壩區(qū)淤積量多424萬m3。

圖7為方案1塔前60 m處淤積高程與出庫流量變形過程。7月1日—8月30日只有發(fā)電洞過流，塔前60 m處持續(xù)淤積抬高，最大值達到181.39 m;9月1日—9月30日除了發(fā)電洞開啟，排沙洞也有開啟，塔前60 m處持續(xù)降低。

圖8為方案2塔前60 m處淤積高程與出庫流量變化過程，7月1日—7月10日發(fā)電洞過流80%，排沙洞過流20%，在此期間來沙較少，塔前60 m處略有淤積抬升，7月11日—8月31日只有發(fā)電洞過流，塔前60 m處持續(xù)淤積抬高，最大值達到181.35 m，9月1日—9月30日除了發(fā)電洞開啟，排沙洞也有開啟，塔前60 m處持續(xù)降低。

方案1與方案2不同之處在于，方案1出庫流量小于1 800 m3/s（6臺機組滿發(fā)）時全部通過發(fā)電洞下泄，方案2出庫流量為400～2 250 m3/s時發(fā)電洞過流80%、排沙洞分流20%。當(dāng)出庫流量大于2 250 m3/s時，兩方案的孔洞運用方式相同。方案2排沙洞開啟的時間略長，有利于排沙。模擬時段內(nèi)底孔前淤積高程均沒有達到允許淤積高程187 m。

3.3.2?排沙比

2003典型7月1日—8月31日進入壩區(qū)的沙量為1.44億t，方案1計算的出庫沙量為0.98億t、排沙比為68.1%，方案2計算的出庫沙量為1.02億t、排沙比為71.0%，方案2比方案1多排沙0.04億t。7月1日—9月30日進入壩區(qū)的沙量為2.31億t，該時段內(nèi)方案1出庫沙量為1.84億t、排沙比為79.9%，方案2出庫沙量為1.89億t、排沙比為81.9%，方案2比方案1多排沙0.05億t，見表5。

3.3.3?發(fā)電效益

水輪機發(fā)電量可表示為

N=9.81ηQH（14）

式中：Q為通過水電站水輪機的流量;H為水電站的凈水頭;η為水電站效率。

經(jīng)計算，2003典型方案1總發(fā)電量為11.57億kW·h，方案2總發(fā)電量為10.65億kW·h，方案1比方案2多發(fā)電0.92億kW·h，見表6。

3.3.4?泄水孔洞調(diào)度方案推薦

綜上，雖然方案1比方案2進水塔前平均淤積厚度多0.91 m，但是兩方案在模擬時段內(nèi)底孔前淤積高程均未達到允許值187 m，都能保證進水塔群前實現(xiàn)防淤堵的要求。在此前提下，方案2底孔開啟次數(shù)多，使得方案2發(fā)電效益不如方案1，為充分發(fā)揮水庫的發(fā)電效益，推薦選用泄水孔洞調(diào)度方案1。

4?結(jié)?語

將建立的立面二維水沙數(shù)學(xué)模型應(yīng)用于小浪底水利樞紐進水塔群的防淤堵研究，在典型設(shè)計水沙條件下，計算和比較了兩種泄水孔洞調(diào)度方案下河床縱剖面的變化及其對排沙比和發(fā)電效益的影響，得出以下結(jié)論。

（1）在模擬時段內(nèi)，方案1和方案2的底孔前淤積高程均未達到允許值187 m。

（2）就平均淤積厚度而言，方案1比方案2在進水塔前多淤積0.91 m;就排沙比而言，在同時段方案2比方案1稍高，但方案2底孔開啟次數(shù)多，發(fā)電效益不如方案1。

（3）在滿足近期進水塔群防淤堵要求的前提下，應(yīng)盡可能發(fā)揮小浪底水庫的發(fā)電效益，泄水孔洞調(diào)度方案推薦選用方案1，即出庫流量小于發(fā)電洞泄流量時優(yōu)先啟用發(fā)電洞泄流，出庫流量超過發(fā)電洞泄流量時超出部分盡量通過排沙洞、明流洞、孔板洞泄流，減少庫區(qū)淤積。

參考文獻：

[1]?王二平，張欣，孫東坡，等.小浪底水利樞紐防泥沙淤堵試驗研究[J].華北水利水電大學(xué)學(xué)報（自然科學(xué)版），2015，36（6）：6-9.

[2]?林秀山.黃河小浪底水利樞紐文集[M].鄭州：黃河水利出版社，2001：1-6.

[3]?邢晨雄.小浪底渾水動床模型驗證及進水塔前沖刷漏斗試驗研究[D].鄭州：華北水利水電大學(xué)，2016：17-22.

[4]?王銀海.高水位水庫淤積泥沙水力排淤模擬研究[D].鄭州：華北水利水電大學(xué)，2018：9-12.

[5]?張欣.小浪底樞紐進水塔前允許淤沙高程試驗研究[D].鄭州：華北水利水電大學(xué)，2016：9-14.

[6]?何勝男.小浪底水利樞紐進水塔群高灘邊岸坍塌特性試驗研究[D].鄭州：華北水利水電大學(xué)，2016：9-14.

[7]?劉子龍，王船海.復(fù)雜邊界準三維水流模擬[J].成都科技大學(xué)學(xué)報，1996，28（1）：12-20.

[8]?李褆來，竇希萍，黃晉鵬.長江口邊界擬合坐標(biāo)的三維潮流數(shù)學(xué)模型[J].水利水運科學(xué)研究，2000（3）：1-6.

[9]?周宜林，道上正規(guī)，檜谷治.非淹沒丁壩附近三維水流運動特性的研究[J].水利學(xué)報，2004，35（8）：46-53.

[10]?余明輝，吳騰，楊國錄.剖面二維水沙數(shù)學(xué)模型及其初步應(yīng)用[J].水力發(fā)電學(xué)報，2006，25（4）：66-69.

[11]?余明輝，楊國錄.平面二維非均勻沙數(shù)值模擬方法[J].水利學(xué)報，2000，31（5）：65-69.

【責(zé)任編輯?張?帥】