基于一維泥沙數學模型的水庫淤積與回水耦合研究

白瑞春

(溫州廣匯水利建設有限公司，浙江 溫州 325011)

0 引言

水庫在我國抗洪抗旱過程中發(fā)揮了巨大作用，全國地表供水量的1/3都由水庫提供.庫容是水庫安全運行的重要保障，而水庫泥沙淤積會導致庫容損失，降低水庫利用效率［1－3］.在各類水庫淤積危害中，庫尾泥沙淤積最為嚴重，由于此處遠離壩前，因此難以對其進行有效沖刷.水庫泥沙淤積還會引起回水抬高，進而引起淹沒、內澇、土地鹽堿化、洪水水位增加等災害，因此研究水庫淤積與回水耦合關系意義重大.研究表明，水庫泥沙運動屬于非均勻不平衡運動［4］.前蘇聯(lián)學者在20世紀30年代首先建立了一維不平衡泥沙輸運方程，為后人的研究奠定了基礎［5］.隨后出現了水沙運動的二維數學模型，通過數值模擬得出了粒徑大于0.1 mm的泥沙會在水庫回水末端沉積［6－7］.基于國內外研究現狀，一維、二維數學模型都已經得到廣泛應用.一維模型適用于長河段河床變形研究;二維模型適用于泥沙運動沿水深、河寬變化不大的流動.結合水動力學、水庫泥沙淤積理論，建立了一維非平衡泥沙輸運數學模型并對其進行數值研究，得出水面線和淤積分布情況.并研究了淤積長度與回水長度、淤積厚度與回水抬升的耦合關系.

1 泥沙運動數學模型的建立

與水沙充分摻混的模型試驗相比，數值模擬不受模型尺寸限制，可以通過改變初值條件、試驗參數等計算多種方案下的泥沙輸運規(guī)律.泥沙運動數學模型主要基于水沙運動規(guī)律，同時將數學模型作為研究手段.常用的數值模擬方法有一維不平衡泥沙模型和二維泥沙模型，本文主要通過一維非平衡泥沙輸運數學模型對水庫泥沙淤積與回水耦合關系進行研究.

依據圣維南定理，懸疑質運動采用一維不平衡輸運方程進行計算［8］.懸疑質流動遵循連續(xù)性方程:

其中:A—河道斷面面積，m2;

Q—河水流量，m3/s;

x—距離，m;

t— 時間，s.

流動還滿足動量方程:

其中:Z—水位，m;g—重力加速度，m3/s2;

n—渠道粗糙系數，無量綱;

R—渠道的水力半徑，m;

qsr—入口流量，m3/s;

ub—來流速度在x方向的分量，m/s.

非平衡泥沙輸運方程為:

其中:β—泥沙濃度系數，無量綱;

Sk—斷面平均泥沙濃度，kg/m3;

ωk—分組泥沙的沉降速度，m3/s;

qsk—出口流量，m3/s;

Ssu，k— 入口平均泥沙濃度，kg/m3;

B—河道斷面寬度，m;

Φk—河道水流攜沙能力，kg/m3.

水庫河床變形方程為:

其中:γ0— 容重，kg/m3;

由此得出第k組泥沙的沖刷函數為:

其中:Pbk—沙床級配函數;

φk—混合沙的挾沙能力.

2 泥沙運動模型的求解

求解一維泥沙數學模型之前，首先需要確定其初始參數，主要參數有:初始粗糙率、泥沙沉降速率、挾沙能力系數;淤積物容重等.經過對試驗對象進行分析，其底坡為J=7.5‰，因此取0～23 m的初始粗糙率為0.017，取23～60 m的初始粗糙率為0.013.假設顆粒直徑為D，泥沙沉降速率按照經驗公式計算得到:

泥沙容重 γs取2650 kg/m3;干容重采 γ取1325 kg/m3.選取一組典型工況，運用Preissmann四點隱式算法和有限差分法對一維泥沙模型進行求解，計算出水面線及淤積分布［6－7］.流量 Q 為0.02 m3/s，尾水水深 H 為 0.44 m，含沙量 S 取1 kg/m3，計算時間長度取7 h.7 h后的水位淤積地形情況(見圖1).

圖1 7 h后的水位淤積地形情況

由圖1可以看出:經過7 h后，在沿程距離10～40 m出已經出現了明顯的泥沙淤積.從試驗結果來看，數值模擬結果與實測結果基本吻合.進一步分析，累計淤積量與沿程距離之間的關系，發(fā)現數值模擬得到的累積淤積量與測得的累積淤積量的偏差在10%以內.7 h后的累積淤積量與沿程距離之間的關系(見圖2).

由圖2可以看出，隨著沿程距離的增加，累積淤積量先平緩上升，隨后快速上升.當沿程距離達到40 m后，上升速度變得十分緩慢.實測數據與數值模擬得到的數據基本吻合，可見本文提出的一維泥沙模型可用于水庫淤積與回水耦合研究.

圖2 累積淤積量與沿程距離的關系

3 水庫淤積與回水耦合研究

由數學模型可以看出，影響淤積的因素很多，但回水抬升是導致水庫淤積的主要原因.在河流上游修筑水庫，庫區(qū)水位抬升，導致水流速度降低，挾沙能力下降，大量的泥沙沉降到水庫中，導致落淤區(qū)域出現回水抬升現象并促進泥沙的進一步淤積［9］.形成泥沙淤積—水位抬升—淤積加劇—回水抬升—再淤積的惡性循環(huán).研究表明，水庫淤積量的大小取決于相對坡比降Jc/J0和相對水深h/H.

其中Jc為淤積平衡比降;J0為原河床比降;H為淤積前的壩前深度;h淤積平衡后的壩前深度.

3.1 淤積長度與回水長度

若底坡為 J=7.5‰，Q 為 0.025 m3/s，尾水水深H為0.45 m，含沙量S取2 kg/m3，沙粒直徑取0.1 mm，計算中的零點為水庫大壩位置.首先計算淤積末端與回水末端隨時間的變化(見圖3).其中圖3上邊的曲線為回水末端，下邊的曲線為淤積末端.

由圖3可知，淤積末端與回水末端隨時間的變化規(guī)律基本相同，但是淤積起始點較回水起始點更靠近水庫大壩.在淤積初步平衡以前，淤積回水末端向上游增長的更快些;但淤積初步平衡后，兩條曲線的增加速度都變得十分緩慢.由于淤積末端和回水末端的增加速度基本一致，因此對其試驗數據做線性擬合，得到計算式:

其中:xwe—回水末端到壩前的距離，m;

xde—淤積末端到壩前的距離，m.

泥沙堆積形成的三角洲洲面淤積長度與回水長度間也存在耦合關系.基于泥沙數學模型對其進行研究，得出泥沙三角洲的情況(見圖4).其中圖4上邊的曲線為回水長度，下邊的曲線為洲面淤積長度.

圖3 淤積末端與回水末端隨時間的變化

由圖4可知，三角洲洲面淤積長度的增長趨勢與回水長度基本一致，兩條曲線呈近似平行關系.回水長度值總是大于洲面長度，在淤積初步平衡之前，兩條曲線的增長速度很快，增加趨勢為近似線性增加;到達初步平衡后，二者增長速度逐漸降低.對二者的試驗數據進行線性擬合，得出本工況回水長度的經驗計算式:

其中:lw—回水長度，m;

lde—三角洲洲面淤積長度，m.

3.2 淤積厚度與回水抬升

在研究淤積長度與回水長度耦合的基礎上，本節(jié)將進一步研究淤積厚度與回水抬升的關系.研究之前，首先定義無量綱距離x/L0，x為到壩前的距離，L0為靜回水長度，本文試驗中的靜回水長度取47.1 m.表1給出了距壩距離、無量綱距離、沿程距離對照表.

表1 沿程、無量綱距離對照表

采用無量綱距離計算各個斷面淤積厚度與水位抬升高度之間的關系，將計算數據繪制成曲線(見圖5).

計算結果顯示，淤積厚度與水位抬升高度的關系十分復雜.圖中:直線adc代表該斷面淤積厚度與水位抬升高度的增加完全相同;曲線bde為初步淤積平衡分界線;直線aef表示淤積三角洲洲頭位置到達各個斷面的時間.由直線adc和曲線bde將圖5劃分為 A、B、C1、C2四個區(qū)域.

區(qū)域A為各個起始淤積點下游斷面未達到初步平衡時的淤積情況.由圖5中可知，各點未達到淤積平衡時，淤積厚度與回水抬升曲線基本與直線adc平行，說明其斜率近似為1.即A區(qū)域內各點淤積厚度對回水的貢獻基本相同.

圖5 淤積厚度、水位抬升高度關系曲線

區(qū)域B為各起始淤積點下游斷面達到初步平衡后的淤積情況.從圖中可以看出，淤積初步平衡后，各點淤積厚度與回水抬升曲線的斜率自左向右依次降低，直到淤積厚度最大處的斜率接近于0.究其原因，此處的壩前深度遠大于淤積厚度，泥沙淤積對水位的影響可以忽略.

區(qū)域C1為各起始淤積點上游斷面達到初步平衡前的淤積情況.區(qū)域C1內的淤積規(guī)律與區(qū)域A相似，此區(qū)域各點淤積厚度對回水的貢獻基本相同，關系直線的斜率接近1.

區(qū)域C2為各起始淤積點上游斷面達到初步平衡后的淤積情況.區(qū)域C2內的淤積規(guī)律與區(qū)域B相似，但是上游斷面淤積厚度對水位抬升的影響較大.

4 結語

庫容是水庫安全運行的重要保障，而水庫泥沙淤積會導致庫容損失，降低水庫利用效率.水庫泥沙淤積還會引起回水抬高，進而引起淹沒、內澇、土地鹽堿化、洪水水位增加等災害，

結合水動力學、水庫泥沙淤積理論，建立了一維非平衡泥沙輸運數學模型并對其進行數值研究，得出水面線和淤積分布情況.并研究了淤積長度與回水長度、淤積厚度與回水抬升的耦合關系.研究表明:

(1)隨著沿程距離的增加，累積淤積量先平緩上升，隨后快速上升.當沿程距離達到40 m后，上升速度變得十分緩慢.實測數據與數值模擬得到的數據基本吻合，可見一維泥沙模型可用于水庫淤積與回水耦合研究.

(2)淤積末端與回水末端隨時間的變化規(guī)律基本相同，但是淤積起始點較回水起始點更靠近水庫大壩.三角洲洲面淤積長度的增長趨勢與回水長度基本一致，兩條曲線呈近似平行關系.

(3)淤積厚度與水位抬升高度的關系十分復雜，本文將其分為4個區(qū)域，更清晰的討論了淤積厚度與水位抬升在不同階段的關系.

［1］盧金友，黃 悅.三峽水庫淤積計算預測與原型實測結果比較分析［J］.長江科學院院報，2013(12):1 －6.

［2］曹慧群，李青云，黃 茁，等.我國水庫淤積防治方法及效果綜述［J］.水力發(fā)電學報，2013(6):183－189.

［3］歐應鈞，封光寅，趙學峰.丹江口水庫泥沙調度方式探討［J］.人民長江，2014(2):82 －85.

［4］張小帥，張耀哲.基于MATLAB的非均勻懸移質不平衡輸沙計算［J］.人民黃河，2014(11):26－29.

［5］李繼偉.梯級水庫群多目標優(yōu)化調度與決策方法研究［D］.保定:華北電力大學，2014.

［6］鄧春艷，吳衛(wèi)民，夏軍強.巖灘水庫泥沙淤積計算數學模型［J］.武漢大學學報:工學版，2013(1):19 －25.

［7］謝金明，吳保生，劉孝盈.水庫泥沙淤積管理綜述［J］.泥沙研究，2013(3):71－80.

［8］謝 龍.三峽變動回水區(qū)末端段復合水動力條件分析及對泥沙輸移的影響［D］.重慶:重慶交通大學，2013.

［9］盧 克.河口懸浮泥沙濃度遙感反演模型的探討［J］.浙江水利水電專科學校學報，2012，73(3):63－65.