由麗華，歐陽莉，趙 陽，王 強(qiáng)，楊 飛

（1.四川省紫坪鋪開發(fā)有限責(zé)任公司，四川 成都 610091；2.黃河水利委員會(huì) 黃河水利科學(xué)研究院，河南 鄭州 450003；3.河南省湖庫功能恢復(fù)與維持工程技術(shù)研究中心，河南 鄭州 450003）

河道內(nèi)挾帶大量細(xì)沙的水流進(jìn)入水庫后，含沙水流因與庫內(nèi)清水存在密度差而潛入水庫底部，形成水庫異重流。水庫異重流對(duì)于水庫的淤積總量以及淤積形態(tài)影響深遠(yuǎn)［1］，因此科學(xué)把握水庫異重流運(yùn)動(dòng)過程可有效判別水庫的淤積形態(tài)，同時(shí)可指導(dǎo)多沙水庫泥沙管理。目前，對(duì)于水庫異重流的研究主要分為理論研究、現(xiàn)場(chǎng)觀測(cè)以及數(shù)學(xué)模型試驗(yàn)等方面。基礎(chǔ)理論研究主要集中在異重流形成條件、潛入點(diǎn)位置及持續(xù)運(yùn)行條件3個(gè)方面。不少學(xué)者通過水槽試驗(yàn)對(duì)異重流的潛入過程進(jìn)行研究，代表學(xué)者包括范家驊［2］、曹如軒等［3］、焦恩澤［4］。Xia等［5］為了研究異重流潛入點(diǎn)判別公式，重點(diǎn)考慮了清水弗勞德數(shù)以及含沙量影響，并通過前人的水槽試驗(yàn)結(jié)果擬合了新的潛入點(diǎn)判別公式。范家驊等［6－7］全面回顧了異重流潛入現(xiàn)象的各種成果，并對(duì)不同進(jìn)口水沙條件下異重流潛入點(diǎn)的水沙特性、潛入點(diǎn)與水沙因子的關(guān)系進(jìn)行了探討?，F(xiàn)階段對(duì)異重流現(xiàn)場(chǎng)觀測(cè)研究多圍繞小浪底水庫開展，自2002年首次開展調(diào)水調(diào)沙試驗(yàn)后［8］，黃河水利委員會(huì)多次開展黃河調(diào)水調(diào)沙運(yùn)行。經(jīng)過多年的調(diào)水調(diào)沙，黃河下游河道主槽不斷萎縮的態(tài)勢(shì)得到了初步遏制，在防洪減淤方面取得了重要進(jìn)展［9］。此后，也有不少學(xué)者針對(duì)黃河調(diào)水調(diào)沙期間小浪底水庫異重流運(yùn)動(dòng)特性進(jìn)行研究，包括潛入點(diǎn)及流速分布［10］、壩前渾液面沉降［11］、壩前異重流厚度［12］以及排沙效率［13］等。由于異重流測(cè)量需要把握時(shí)機(jī)，且現(xiàn)階段其他水庫較少有現(xiàn)場(chǎng)觀測(cè)的資料，因此構(gòu)建數(shù)學(xué)模型研究水庫異重流運(yùn)動(dòng)過程也是目前較為主流的研究手段之一。熊正偉［14］構(gòu)建了同時(shí)適用于異重流及溯源沖刷的水沙輸移與河床沖淤耦合的一維數(shù)學(xué)模型，并對(duì)1975年三門峽水庫溯源沖刷過程進(jìn)行了模擬，分析了溯源沖刷發(fā)展特點(diǎn)及規(guī)律。常云華等［15］提出采用一維非飽和不平衡輸沙異重流數(shù)學(xué)模型進(jìn)行水庫泥沙淤積計(jì)算，在黃河潼關(guān)至三門峽庫區(qū)段進(jìn)行計(jì)算的結(jié)果良好。王增輝［1］綜合考慮干支流倒灌因素，對(duì)水庫溯源沖刷以及異重流過程進(jìn)行模擬，同時(shí)改進(jìn)淺水二維水沙耦合模型對(duì)不同邊界條件下的水庫支流口溯源沖刷進(jìn)行模擬研究。筆者針對(duì)近年來紫坪鋪水庫來沙偏多、水庫淤積日漸嚴(yán)重的現(xiàn)狀，開展異重流與明流耦合水沙數(shù)學(xué)模型構(gòu)建工作，并利用2019年實(shí)測(cè)水沙數(shù)據(jù)進(jìn)行驗(yàn)證，以期為紫坪鋪水庫開展水庫清淤與泥沙調(diào)度提供指導(dǎo)。

1 研究區(qū)域與數(shù)據(jù)

紫坪鋪水利樞紐工程位于四川省岷江上游映秀鎮(zhèn)至都江堰河段，是一座以灌溉和供水為主，兼有防洪、發(fā)電、環(huán)境保護(hù)和旅游等綜合效益的一等大（I）型水利水電工程，也是都江堰灌區(qū)和成都市的重要水源調(diào)節(jié)工程。壩址以上流域面積2.26萬km2，多年平均徑流量148億m3，占岷江上游徑流總量的97%。水庫正常蓄水位877.0 m，相應(yīng)庫容9.98億m3；校核洪水位883.1 m，相應(yīng)庫容11.12億m3；死水位817.0 m，相應(yīng)庫容2.24億m3；水庫調(diào)節(jié)庫容7.74億m3。水庫設(shè)計(jì)多年平均懸移質(zhì)沙量792萬t，含沙量0.572 kg／m3。多年汛期（6—10月）平均沙量778.9萬t，占全年沙量的98.3%；5—9月沙量為760.7萬t，占全年沙量的96%；7—9月沙量為533.43萬t，占全年沙量的67.4%。多年平均入庫推移質(zhì)70萬t，最大粒徑為500 mm，中值粒徑為113 mm。紫坪鋪水庫入庫控制站為映秀水文站，出庫控制站為紫坪鋪水文站。

紫坪鋪水庫2019年來水偏多，形成幾次較大的洪峰，最大洪峰為1 792 m3／s。受降雨影響，上游發(fā)生多起泥石流以及滑坡現(xiàn)象，造成罕見的入庫沙峰，最大沙峰為63.4 kg／m3。圖1為紫坪鋪水庫2019年汛期進(jìn)出庫流量與水位變化過程，由圖1可知，2019年汛期（6—10月）水量總體偏大，形成幾次不小的洪峰，水位按照調(diào)度規(guī)程，圍繞850 m調(diào)度。然而，如圖2所示，入庫泥沙只出現(xiàn)了一次較大的峰值（8月20日），出庫含沙量出現(xiàn)了歷史罕見的63.4 kg／m3。由于2019年有異重流到達(dá)壩前，因此選用2019年的數(shù)據(jù)進(jìn)行模型驗(yàn)證。

圖1 紫坪鋪水庫2019年汛期進(jìn)出庫流量與水位變化過程

圖2 紫坪鋪水庫2019年汛期進(jìn)出庫含沙量變化過程

2 一維水沙數(shù)學(xué)模型構(gòu)建

基于水流連續(xù)性方程、水流動(dòng)量方程、泥沙連續(xù)方程、河床變形方程以及異重流運(yùn)動(dòng)方程構(gòu)建紫坪鋪水庫一維水沙數(shù)學(xué)模型。天然水庫的邊界條件十分復(fù)雜，為了便于求解，本研究對(duì)部分條件做如下簡(jiǎn)化：①假定可以將最小計(jì)算時(shí)段內(nèi)的非恒定流作為恒定流處理，具體是將進(jìn)口斷面的實(shí)際流量過程線改為若干個(gè)不同流量級(jí)組成的梯級(jí)過程線進(jìn)行計(jì)算，對(duì)于每一個(gè)梯級(jí)來說，流量為常數(shù)，水流為恒定流；②假定河床沖淤過程中每一個(gè)短時(shí)段內(nèi)河床變形對(duì)水流條件影響不大；③在最小計(jì)算時(shí)段內(nèi)水體中含沙量恒定。天然庫區(qū)的泥沙是非均勻沙，因此在實(shí)際的模擬計(jì)算中，需要將泥沙顆粒級(jí)配按照粒徑大小進(jìn)行分級(jí)處理，最終按照分級(jí)級(jí)配計(jì)算水流的挾沙能力。其中，水流泥沙基本運(yùn)動(dòng)方程見文獻(xiàn)［15］，異重流模塊主要分為計(jì)算潛入點(diǎn)的水深、異重流均勻流水深、判斷異重流發(fā)生條件以及異重流的運(yùn)動(dòng)方程，具體計(jì)算過程如下。

（1）計(jì)算潛入點(diǎn)的水深h0。公式為

式中：h0為潛入點(diǎn)水深；Q為流量；ηg為重力修正系數(shù)；g為重力加速度；B為河寬；S為含沙量；γs為泥沙容重，取2 650 kg／m3；γ0為清水容重，取1 000 kg／m3。

（2）計(jì)算異重流均勻流水深hn。公式為

式中：hn為均勻流水深；λ為異重流阻力系數(shù)，取λ＝0.025～0.03；J為底坡（實(shí)際為能坡）。

（3）判斷水深h。當(dāng)異重流潛入且能以均勻流形式運(yùn)動(dòng)時(shí)，其潛入點(diǎn)水深h＞max（h0，hn），若滿足，則可以認(rèn)定異重流潛入，并能以均勻流方式運(yùn)行，不上升到水面（式中：A為過水面積；i－1和i分別為前后兩個(gè)計(jì)算時(shí)間點(diǎn)）。

根據(jù)上述方程以及水流泥沙基本運(yùn)動(dòng)方程，將紫坪鋪水庫一維水沙數(shù)學(xué)模型的計(jì)算步驟總結(jié)如下：①輸入邊界條件、初始條件及有關(guān)參數(shù)；②計(jì)算各水力要素；③計(jì)算不同斷面的平均含沙量；④計(jì)算子斷面含沙量；⑤計(jì)算子斷面的沖淤面積；⑥修正子斷面的河床組成；⑦按照相鄰子斷面的沖淤厚度以及過水?dāng)嗝娴膶挾人嫉臋?quán)重，修正節(jié)點(diǎn)上的河底高程；⑧輸出計(jì)算結(jié)果，進(jìn)行下一斷面的計(jì)算。

3 模型結(jié)果及分析

3.1 出庫流量模擬與分析

由于紫坪鋪水庫于2018年才開始建設(shè)映秀水文站，因此在此之前并無實(shí)測(cè)的入庫含沙量資料。本文選取2019年的水沙、地形資料對(duì)一維異重流水沙數(shù)學(xué)模型進(jìn)行驗(yàn)證。具體的參數(shù)設(shè)置如下：①地形數(shù)據(jù)，選擇大斷面數(shù)據(jù)（不包括加測(cè)斷面數(shù)據(jù)），斷面地形經(jīng)過內(nèi)插后，精度達(dá)到1 m；②大斷面的距壩里程，按照河道深泓線進(jìn)行量測(cè)；③不可沖刷斷面，由于紫坪鋪水庫屬于逐漸淤積型水庫，因此不可沖刷地形設(shè)置為2008年的初始地形，認(rèn)為無論如何沖刷，地形均比2008年初始地形要高；④進(jìn)口水沙數(shù)據(jù)，以2019年的進(jìn)出庫水沙為例，進(jìn)口控制站使用映秀水文站的數(shù)據(jù)，出庫站使用紫坪鋪水文站的數(shù)據(jù)；⑤壩前水位，按照紫坪鋪水庫2019年壩前水位進(jìn)行計(jì)算；⑥支流庫容，采用紫坪鋪水庫2019年實(shí)測(cè)支流庫容數(shù)據(jù)；⑦床沙級(jí)配，采用2020年實(shí)測(cè)河床級(jí)配；⑧懸沙粒徑，以取樣分析結(jié)果為主。

圖3為紫坪鋪水庫實(shí)測(cè)出庫流量與計(jì)算出庫流量對(duì)比，圖中黑色線條表示計(jì)算出庫流量過程，紅色線條表示實(shí)測(cè)出庫流量過程。由圖3可知，兩條曲線的吻合度較高。流量計(jì)算結(jié)果與實(shí)測(cè)數(shù)據(jù)非常吻合，無論是峰值還是過程線都達(dá)到了較高的精度，含沙量計(jì)算結(jié)果與紫坪鋪水文站實(shí)測(cè)結(jié)果也較為接近，這也反映出模型的參數(shù)選擇基本上可以作為紫坪鋪水庫水沙運(yùn)動(dòng)情景計(jì)算的參數(shù)。

圖3 紫坪鋪水庫實(shí)測(cè)與計(jì)算出庫流量變化過程對(duì)比

3.2 出庫含沙量模擬與分析

2019年，紫坪鋪水庫遭遇了“8·20”特大山洪泥石流災(zāi)害，紫坪鋪入庫流量大幅上漲，由20日3時(shí)的884 m3／s猛漲到4時(shí)的3 457 m3／s，洪水過程形成迅速。水庫水位也持續(xù)上漲，從3時(shí)的845.19 m上漲到8時(shí)的848.13 m。為應(yīng)對(duì)洪水，水庫管理部門開始通過調(diào)度方式對(duì)洪水進(jìn)行攔蓄，將最大洪峰3 457 m3／s削減至620 m3／s出庫，削峰率為82%。至20日16時(shí)，紫坪鋪水庫共攔蓄洪水6 266萬m3，有效緩解了下游河道防洪壓力，保護(hù)了下游城市安全。在本次調(diào)度過程中，水庫管理部門于20日8時(shí)50分在水位848.75 m時(shí)開啟沖沙放空洞泄洪，14時(shí)水位迫近汛限水位850 m時(shí)開啟1＃泄洪排沙洞進(jìn)行泄洪，并控制最大出庫流量不超過1 800 m3／s。上游暴雨山洪泥石流引發(fā)高濁度洪水進(jìn)入水庫，導(dǎo)致紫坪鋪泄水濁度達(dá)到30 000 NTU，16時(shí)沖沙放空洞全關(guān)。

為驗(yàn)證2019年紫坪鋪水庫汛期含沙量變化過程，本文利用一維異重流水沙模型計(jì)算出庫含沙量，并與出庫站紫坪鋪水文站含沙量變化過程進(jìn)行對(duì)比。圖4為計(jì)算與實(shí)測(cè)出庫含沙量變化過程對(duì)比，整體上看，兩條過程線吻合程度較高。實(shí)測(cè)最大出庫含沙量為7.73 kg／m3，計(jì)算最大出庫含沙量為7.673 kg／m3，吻合度較高，其他時(shí)段的計(jì)算出庫含沙量略大于實(shí)測(cè)出庫含沙量，究其原因，主要是低含沙量的泥沙即使以異重流的形式到達(dá)壩前，若排沙設(shè)施不打開，出庫含沙量依舊很低，其出庫過程與水庫調(diào)度過程相關(guān)。

圖4 紫坪鋪水庫實(shí)測(cè)與計(jì)算出庫含沙量變化過程對(duì)比

3.3 庫區(qū)深泓線模擬與分析

深泓線即沿河流方向最大水深處的連線，沿此線的剖面為河流的縱剖面。深泓線的高程基本可以反映河床縱向變化情況。圖5為模型模擬2019年汛后與2008年汛后實(shí)測(cè)、2019年汛前實(shí)測(cè)以及2020年汛后實(shí)測(cè)深泓線對(duì)比。由圖5可知，2019年汛期結(jié)束后的深泓線基本在2019年汛前與2020年汛后之間，因此本文構(gòu)建的一維異重流水沙數(shù)學(xué)模型的精度相對(duì)較高，既可以對(duì)流量過程進(jìn)行精確模擬，也可以對(duì)輸沙過程進(jìn)行精確模擬。

圖5 紫坪鋪水庫實(shí)測(cè)與計(jì)算深泓線變化過程對(duì)比

3.4 庫區(qū)淤積分布模擬與分析

利用一維異重流水沙數(shù)學(xué)模型可以模擬水庫不同位置的淤積量，圖6為紫坪鋪水庫2019年汛期不同位置單位距離實(shí)測(cè)與計(jì)算淤積量結(jié)果對(duì)比，黑色實(shí)線為模型計(jì)算結(jié)果，紅色實(shí)線為實(shí)測(cè)結(jié)果，兩條曲線的趨勢(shì)大致相同，淤積的位置相對(duì)準(zhǔn)確。然而，紅色實(shí)線較黑色實(shí)線總體偏高，也就是說，實(shí)測(cè)的淤積量比模擬計(jì)算的淤積量大，究其原因，主要是模型計(jì)算結(jié)果僅為2019年汛期的淤積量，實(shí)測(cè)淤積量為2019—2020年兩年汛期的淤積量，由于2020年水沙數(shù)據(jù)未獲取，因此并沒有復(fù)核2020年紫坪鋪水庫的淤積形態(tài)。但是通過實(shí)測(cè)與計(jì)算淤積量的對(duì)比可知，該模型可以計(jì)算庫區(qū)不同位置的淤積量，且具有較好的精度。由圖6可知，2019年模型計(jì)算結(jié)果為汛期淤積2 274萬m3，而2019—2020年實(shí)測(cè)淤積量為2 477萬m3。2020年汛期水量比較大，但是來沙量相對(duì)偏少，可以推測(cè)，2020年全年的淤積量較2019年少很多，因此該模型計(jì)算的結(jié)果較為合理。

圖6 紫坪鋪水庫2019年汛期不同位置單位距離實(shí)測(cè)與計(jì)算淤積量結(jié)果對(duì)比

4 結(jié) 論

（1）基于水流連續(xù)性方程、水流動(dòng)量方程、泥沙連續(xù)方程、河床變形方程以及異重流運(yùn)動(dòng)方程構(gòu)建紫坪鋪一維水沙數(shù)學(xué)模型。依據(jù)實(shí)測(cè)資料對(duì)該模型進(jìn)行了率定與參數(shù)的確定，最終得到了可以計(jì)算異重流演進(jìn)過程的水沙一維數(shù)學(xué)模型。

（2）汛期的流量數(shù)據(jù)驗(yàn)證結(jié)果表明：流量計(jì)算結(jié)果與實(shí)測(cè)數(shù)據(jù)較為吻合，含沙量計(jì)算結(jié)果與紫坪鋪水文站實(shí)測(cè)結(jié)果也較為接近，這反映出本模型的參數(shù)選擇基本上可以作為紫坪鋪水庫水沙運(yùn)動(dòng)情景計(jì)算的參數(shù)。

（3）汛期的含沙量數(shù)據(jù)驗(yàn)證結(jié)果表明：實(shí)測(cè)最大出庫含沙量誤差很小，吻合度較高，其他時(shí)段的計(jì)算出庫含沙量略大于實(shí)測(cè)出庫含沙量。究其原因，主要是因?yàn)榈秃沉康哪嗌臣词挂援愔亓鞯男问降竭_(dá)壩前，若排沙設(shè)施不打開，出庫含沙量依舊很低，其出庫過程與水庫調(diào)度過程相關(guān)。

（4）汛期的地形數(shù)據(jù)驗(yàn)證結(jié)果表明：2019年汛期結(jié)束后的深泓線基本在2019年汛前與2020年汛后之間，水沙數(shù)學(xué)模型的精度相對(duì)較高。2019年模型計(jì)算結(jié)果為汛期淤積2 274萬m3，而2019—2020年實(shí)測(cè)淤積量為2 477萬m3。