江穎 李珍 安靜泊

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

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

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

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

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

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

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

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

水流連續(xù)方程：

ux+vy=0（1）

水流運(yùn)動(dòng)方程：

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

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

懸移質(zhì)運(yùn)動(dòng)方程：

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

河床變形方程：

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

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

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

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

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

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

2.1?方程離散

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

水流運(yùn)動(dòng)方程可離散為

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為時(shí)間步長(zhǎng);n為時(shí)間節(jié)點(diǎn)編號(hào);i為x方向節(jié)點(diǎn)編號(hào);j為y方向節(jié)點(diǎn)編號(hào)。

懸移質(zhì)運(yùn)動(dòng)方程在矩形網(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）右端項(xiàng)的離散代數(shù)關(guān)系式，根據(jù)文獻(xiàn)[11]進(jìn)行求解。

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?封閉模式

泥沙運(yùn)動(dòng)極其復(fù)雜，對(duì)于其中河底挾沙力公式需要采用經(jīng)驗(yàn)公式來(lái)進(jìn)行封閉。立面二維水沙數(shù)學(xué)模型計(jì)算中與河床沖淤變形相關(guān)的因素是河床底部含沙量與挾沙力，通常選一接近河底的點(diǎn)作為參考點(diǎn)，參考點(diǎ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分別為計(jì)算參數(shù)，j1=∫1ξa[（1-ξ）/ξ]Zdξ，j2=∫1ξa[（1-ξ）/ξ]Zln ξdξ，其中Z為懸浮指標(biāo)、ξ為水深;k為卡門常數(shù);z1=ωβkU*=Zβ，其中β為參數(shù)、U*為摩阻流速。

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

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

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

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

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

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

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

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

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

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

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

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

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

3.3.1?河床縱剖面分析

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

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

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

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

3.3.2?排沙比

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

3.3.3?發(fā)電效益

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

N=9.81ηQH（14）

式中：Q為通過(guò)水電站水輪機(jī)的流量;H為水電站的凈水頭;η為水電站效率。

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

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

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

4?結(jié)?語(yǔ)

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

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

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

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

參考文獻(xiàn)：

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

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

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

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

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

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

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

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

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

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

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

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