基于滯后響應(yīng)模型的三門峽水庫(kù)沖淤計(jì)算方法

沈 逸，鄭 珊，吳保生

(1.清華大學(xué) 水沙科學(xué)與水利水電工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100084；2.武漢大學(xué) 水資源與水電工程科學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室，湖北 武漢 430072)

1 研究背景

沖積河流與流域不斷進(jìn)行著物質(zhì)和能量交換。由于河流因素的復(fù)雜性，上游輸入的泥沙一般難以與當(dāng)前河段的水流挾沙力相等，導(dǎo)致沖積河流的河床時(shí)刻處于沖淤變形中。然而，沖積河流的自動(dòng)調(diào)整作用又使得河床演變朝著變形停止的方向發(fā)展，顯示出平衡傾向性[1]。

河床受到外界擾動(dòng)后的自動(dòng)調(diào)整，其初始速度通常較快，但隨著河床不斷趨近于新的平衡狀態(tài)，調(diào)整速率逐漸減緩，并最終趨近于零，如圖1所示。宏觀尺度上的河床變形總是滯后于來(lái)水來(lái)沙條件等外部因素的變化，稱為河床演變的滯后響應(yīng)[2]。對(duì)于圖1所示的滯后響應(yīng)模式，吳保生等[3-4]基于變率公式建立了河床滯后響應(yīng)模型的理論框架，包括通用積分、單步解析、多步遞推三種計(jì)算模式，適用于模擬不同條件下的河床演變過(guò)程。其中，滯后響應(yīng)模型的單步解析公式可表示為：

y=(1-e-βt)ye+e-βty0

(1)

式中：y為特征變量當(dāng)前值；ye為特征變量平衡值；y0為特征變量初始值；β為調(diào)整速率參數(shù)；t為調(diào)整時(shí)間。

圖1 滯后響應(yīng)模式及特征變量準(zhǔn)平衡時(shí)間示意

特征變量平衡值ye和調(diào)整速率參數(shù)β是滯后響應(yīng)模型的兩個(gè)關(guān)鍵參數(shù)，平衡值ye表示河道受擾動(dòng)后，特征變量調(diào)整的目標(biāo)值；調(diào)整速率參數(shù)β則表示特征變量調(diào)整的快慢，β越大表示特征變量調(diào)整速率越快。平衡值ye往往根據(jù)具體研究對(duì)象的特點(diǎn)給出計(jì)算方法。例如，張為等[5]基于水流挾沙力公式推導(dǎo)得到挾沙強(qiáng)度的計(jì)算式，根據(jù)挾沙強(qiáng)度與汛期平均含沙量的比值反映河道的沖淤過(guò)程，得到三峽水庫(kù)泥沙淤積量平衡值的計(jì)算公式。呂宜衛(wèi)等[6]在推導(dǎo)下荊江河段單位河長(zhǎng)累計(jì)沖刷量平衡值計(jì)算公式時(shí)，在考慮挾沙強(qiáng)度和汛期來(lái)沙量影響的基礎(chǔ)上，引入了無(wú)量綱水面比降來(lái)減弱洞庭湖頂托效應(yīng)對(duì)模型計(jì)算的影響。鄭珊[7]將擾動(dòng)后河道淤積體的平衡縱剖面形態(tài)概化為三角形或梯形，為累計(jì)淤積量平衡值的計(jì)算提供了一種概化方法。調(diào)整速率參數(shù)β多取為常數(shù)，僅少數(shù)研究在模擬黃河下游平灘流量時(shí)，把調(diào)整速率參數(shù)β表示為汛期平均流量的冪函數(shù)[8]。

此外，關(guān)于河床調(diào)整的滯后響應(yīng)時(shí)間，以往研究主要基于資料相關(guān)分析，缺乏嚴(yán)格的數(shù)學(xué)或理論推導(dǎo)。例如，張?jiān)h等[9]發(fā)現(xiàn)潼關(guān)高程與6年滑動(dòng)平均年水量具有較好的相關(guān)關(guān)系。廖治棋等[10]發(fā)現(xiàn)監(jiān)利站平灘面積受到前期3年內(nèi)來(lái)水來(lái)沙條件的影響較大。唐小婭等[11]通過(guò)計(jì)算發(fā)現(xiàn)三峽水庫(kù)汛期泥沙累計(jì)淤積量與5年線性疊加壩前水位關(guān)系最密切。王彥君等[12]的研究結(jié)果表明，黃河下游主槽斷面形態(tài)受當(dāng)年在內(nèi)的前8年水沙條件的累積影響。

河床演變的滯后響應(yīng)特征在水庫(kù)淤積過(guò)程中體現(xiàn)明顯。三門峽水庫(kù)自1960年9月蓄水運(yùn)用以來(lái)，泥沙淤積一直是制約水庫(kù)防洪和發(fā)電等綜合效益發(fā)揮的關(guān)鍵問(wèn)題。過(guò)去圍繞水庫(kù)淤積和潼關(guān)高程變化與控制開(kāi)展了大量研究[13-15]。吳保生等[16]發(fā)現(xiàn)三門峽庫(kù)區(qū)淤積不僅與當(dāng)年的流量和運(yùn)行條件有關(guān)，還受到前3～4年的流量和運(yùn)行條件的影響。近期吳保生等[17]采用水沙參數(shù)的多年平均值計(jì)算參數(shù)β，難以體現(xiàn)當(dāng)前年河床沖淤狀態(tài)對(duì)調(diào)整速率參數(shù)的影響。

本文將首先基于水庫(kù)沖淤平衡縱剖面幾何關(guān)系及其物理圖景，推導(dǎo)得到潼關(guān)高程和庫(kù)區(qū)累計(jì)淤積量平衡值的計(jì)算模式。然后引入調(diào)整速率參數(shù)β隨河道沖淤及水沙條件發(fā)生變化的計(jì)算方法，采用式(1)表示的滯后響應(yīng)模型單步解析模式，建立庫(kù)區(qū)累計(jì)淤積量及潼關(guān)高程的計(jì)算方法。進(jìn)而計(jì)算并分析河道達(dá)到準(zhǔn)平衡狀態(tài)所需要的時(shí)間，對(duì)比河道實(shí)際完成調(diào)整量與目標(biāo)調(diào)整量的差距，揭示河床演變的滯后響應(yīng)規(guī)律。

2 潼關(guān)高程與累計(jì)淤積量滯后響應(yīng)模型公式推導(dǎo)

2.1 潼關(guān)高程平衡值三門峽水庫(kù)自1973年11月以來(lái)采用蓄清排渾的運(yùn)用方式，非汛期抬高水位運(yùn)用，汛期平水期控制水位305 m發(fā)電，遇洪水時(shí)降低水位泄洪排沙。胡春宏等[18]指出三門峽水庫(kù)蓄清排渾的運(yùn)用方式是一個(gè)根據(jù)來(lái)水來(lái)沙條件和實(shí)際需求不斷調(diào)整的過(guò)程。自蓄清排渾運(yùn)用開(kāi)始到1985年期間，汛期4個(gè)月全部敞泄排沙，后又發(fā)展為當(dāng)潼關(guān)站出現(xiàn)3000 m3/s流量的洪水期時(shí)水庫(kù)敞泄排沙。1986—2000年間，由于汛期大于3000 m3/s流量的天數(shù)減少，敞泄排沙的入庫(kù)流量標(biāo)準(zhǔn)下調(diào)至2500 m3/s。2000年之后，敞泄排沙的入庫(kù)流量標(biāo)準(zhǔn)又進(jìn)一步降至2000～1500 m3/s。

圖2 三門峽水庫(kù)深泓縱剖面示意

圖2顯示1974年汛后河道深泓縱剖面可用直線擬合，擬合直線的斜率可表示為深泓縱剖面比降Js，截距即為壩前斷面深泓高程Z*s。據(jù)此，可以得到圖3所示水庫(kù)沖淤平衡狀態(tài)下的水庫(kù)縱剖面示意圖，潼關(guān)高程平衡值Ze(m)及潼關(guān)站深泓高程Ztgs(m)分別表示為：

Ze=Ztgs+htg

(2)

Ztgs=Z*es+JesL

(3)

式中：htg為潼關(guān)斷面在1000 m3/s流量時(shí)的水深，m；Z*es為水庫(kù)沖淤平衡時(shí)河道深泓縱剖面起始高程，即壩前斷面深泓高程，m；Jes為水庫(kù)沖淤平衡時(shí)河道深泓縱剖面比降；L為潼關(guān)斷面至大壩的河流長(zhǎng)度，取113.5 km。

(4)

(5)

式中：Qtgi為日均入庫(kù)流量；Si為日均入庫(kù)含沙量；Zdi為日均壩前水位；n為加權(quán)計(jì)算的天數(shù)，平年n=365，閏年n=366；ε和σ分別為流量和含沙量的指數(shù)。

根據(jù)以往研究[19-20]，河道深泓縱剖面比降Jes與水沙條件及水庫(kù)運(yùn)用水位相關(guān)，可表示為：

(6)

將式(4)和式(6)代入式(3)，得到Ztgs的表達(dá)式：

(7)

進(jìn)而將式(7)代入式(2)，即可得到潼關(guān)高程平衡值Ze的完整表達(dá)式：

(8)

式(8)可進(jìn)一步簡(jiǎn)化為：

(9)

(10)

圖3 三門峽水庫(kù)沖淤平衡縱剖面概化示意

2.2 累計(jì)淤積量平衡值為了簡(jiǎn)化計(jì)算，將河道斷面形態(tài)概化為寬度為B的矩形，進(jìn)而可以得到由初始河床縱剖面和與當(dāng)前河床縱剖面組成的概化水庫(kù)平衡淤積體?？紤]到淤積量計(jì)算中的概化河底高程與深泓高程存在一定差別，圖3中把累計(jì)淤積量平衡值計(jì)算中涉及的概化河底高程和河床縱剖面的符號(hào)統(tǒng)一用下標(biāo)m表示。

累計(jì)淤積量平衡值計(jì)算中涉及的初始河床邊界參數(shù)，有概化河底縱剖面的初始比降Jm0、潼關(guān)斷面的初始概化河底高程Ztgm0和壩前斷面的初始概化河底高程Zdm0。

圖3中的沖淤平衡時(shí)的概化河底縱剖面比降、概化河底縱剖面起始高程及潼關(guān)斷面概化河底高程分別記為Jem、Z*em、Ztgm。類似式(4)和(7)，Z*em和Ztgm可分別表示為：

(11)

(12)

(13)

按照?qǐng)D3所示水庫(kù)沖淤平衡縱剖面，累計(jì)淤積量平衡值Ve可表示為：

(14)

分別將Z*em和Ztgm的計(jì)算式(11)和(12)代入式(14)有：

(15)

進(jìn)一步化簡(jiǎn)，可得如下累計(jì)淤積量平衡值計(jì)算公式：

(16)

2.3 滯后響應(yīng)模型單步解析模式迭代計(jì)算方法當(dāng)考慮β隨時(shí)段變化時(shí)，應(yīng)用滯后響應(yīng)模型的多步遞推模式計(jì)算較為復(fù)雜[17]，因此，本文采用單步解析模式(式(1))對(duì)潼關(guān)高程和累計(jì)淤積量進(jìn)行計(jì)算。對(duì)于外部擾動(dòng)階梯狀變化的情況，一個(gè)時(shí)段(水文年)的河床調(diào)整結(jié)果，無(wú)論是否已經(jīng)達(dá)到平衡狀態(tài)，都將作為下一個(gè)時(shí)段的初始條件對(duì)河床演變產(chǎn)生影響，并由此使前期的水沙條件對(duì)后期的河床演變產(chǎn)生影響[2]。

因此，在本文模型中，每個(gè)時(shí)段的特征變量平衡值分別根據(jù)對(duì)應(yīng)時(shí)段內(nèi)的年均入庫(kù)流量、年均來(lái)沙系數(shù)、年加權(quán)水位和年均壩前水位確定，將上一時(shí)段特征變量的計(jì)算結(jié)果作為下一時(shí)段的初始條件，采用式(1)逐時(shí)段計(jì)算每個(gè)時(shí)段末特征變量的狀態(tài)值，依次遞推便可以得到經(jīng)過(guò)多個(gè)時(shí)段后的特征變量狀態(tài)值，如圖4。

圖4 單步解析模式迭代計(jì)算流程示意

單步遞推公式的一般形式可表示為：

yi=(1-e-βΔt)ye，i+e-βΔtyi-1

(17)

式中：Δt為時(shí)段長(zhǎng)度；i為時(shí)段序號(hào)。

將潼關(guān)高程和累計(jì)淤積量作為特征變量代入式(17)，通過(guò)逐時(shí)段遞推可以得到潼關(guān)高程和累計(jì)淤積量的逐年計(jì)算值：

Zi=(1-e-βΔt)Ze，i+e-βΔtZi-1

(18)

Vi=(1-e-βΔt)Ve，i+e-βΔtVi-1

(19)

式中：Zi-1和Zi分別為第i-1和i年末的潼關(guān)高程計(jì)算值；Ze，i為第i年的潼關(guān)高程平衡值；Z0為潼關(guān)高程初始值，采用1974年末潼關(guān)高程實(shí)測(cè)值；Vi-1和Vi分別為第i-1和i年末的累計(jì)淤積量計(jì)算值；Ve，i為第i年的累計(jì)淤積量平衡值；V0為累計(jì)淤積量初始值，采用1966年末累計(jì)淤積量實(shí)測(cè)值。在潼關(guān)高程和累計(jì)淤積量的計(jì)算中，Δt=1 a，調(diào)整速率參數(shù)β的單位為a-1。式(18)和式(19)中特征變量平衡值Ze，i和Ve，i分別采用式(9)和式(16)計(jì)算。

調(diào)整速率參數(shù)β需綜合考慮來(lái)水來(lái)沙條件及河道沖淤狀態(tài)的影響。河道調(diào)整速率與水流挾沙力及水體中含沙量大小相關(guān)，當(dāng)水體中床沙質(zhì)含沙量大于挾沙力時(shí)，泥沙將落淤，含沙量越大，河床淤積越顯著；反之，當(dāng)水體中床沙質(zhì)含沙量小于挾沙力時(shí)，沖刷河床，沖刷與水動(dòng)力強(qiáng)弱緊密相關(guān)[21]。此外，以往研究表明河流受擾動(dòng)后的調(diào)整速率往往先快后慢[22]。綜上所述，假設(shè)庫(kù)區(qū)淤積或潼關(guān)高程抬升時(shí)(即平衡值大于初始值)，β為含沙量的指數(shù)衰減函數(shù)；當(dāng)庫(kù)區(qū)沖刷或潼關(guān)高程降低時(shí)(即平衡值小于初始值)，β為流量的指數(shù)衰減函數(shù)：

(20)

3 潼關(guān)高程和累計(jì)淤積量的模擬與驗(yàn)證

3.1 潼關(guān)高程和累計(jì)淤積量計(jì)算結(jié)果根據(jù)1975—2018年實(shí)測(cè)水沙資料和潼關(guān)高程數(shù)據(jù)，通過(guò)多元非線性回歸方法，得到潼關(guān)高程單步遞推公式(式(18))中平衡值Ze和調(diào)整速率參數(shù)βtg計(jì)算公式的參數(shù)值。然后，使用擬合得到的單步遞推公式(式(18))計(jì)算1975—2018年的潼關(guān)高程。

采用多元回歸擬合得到潼關(guān)高程平衡值Ze和調(diào)整速率參數(shù)βtg的計(jì)算公式分別為：

(21)

(22)

(23)

由式(22)可得潼關(guān)高程處于上升狀態(tài)時(shí)調(diào)整速率參數(shù)與年均含沙量正相關(guān)；潼關(guān)高程處于下降狀態(tài)時(shí)調(diào)整速率參數(shù)βtg與年均入庫(kù)流量正相關(guān)，但隨著流量的增大調(diào)整速率參數(shù)增加不明顯。

潼關(guān)高程歷年計(jì)算值與實(shí)測(cè)值的對(duì)比情況見(jiàn)圖5(a)。潼關(guān)高程實(shí)測(cè)值與計(jì)算值在1975—2002年的決定系數(shù)R2為0.80，1980—2002年以及1980—2018年的決定系數(shù)R2均為0.88，在2003—2018年的決定系數(shù)R2為0.75。計(jì)算值與實(shí)測(cè)值在1975—2002年和2003—2018年相對(duì)誤差的平均值分別為0.06%和0.02%，均處于較低水平。由圖5(b)可知，潼關(guān)高程計(jì)算值與實(shí)測(cè)值接近，計(jì)算效果較以往有所提高。其中，1980—2018年和2003—2018年的決定系數(shù)R2分別從之前的最好計(jì)算結(jié)果0.85和0.74[17]提升至0.88和0.75。

根據(jù)1967—2018年實(shí)測(cè)水沙資料和累計(jì)淤積量數(shù)據(jù)，通過(guò)多元回歸方法得到累計(jì)淤積量單步遞推公式(式(19))中關(guān)鍵參數(shù)平衡值Ve和調(diào)整速率參數(shù)βV的計(jì)算公式。然后，使用擬合得到的單步遞推公式(式(19))計(jì)算1967—2018年的累計(jì)淤積量。擬合得到累計(jì)淤積量平衡值Ve和調(diào)整速率參數(shù)βV的計(jì)算公式分別為：

(24)

(25)

(26)

由式(25)可知庫(kù)區(qū)處于淤積狀態(tài)時(shí)調(diào)整速率參數(shù)與年均含沙量正相關(guān)；庫(kù)區(qū)處于沖刷狀態(tài)時(shí)調(diào)整速率參數(shù)與年均入庫(kù)流量正相關(guān)。

圖5 潼關(guān)高程計(jì)算值與實(shí)測(cè)值對(duì)比

累計(jì)淤積量歷年計(jì)算與實(shí)測(cè)值的對(duì)比情況見(jiàn)圖6(a)。計(jì)算值與實(shí)測(cè)值在各個(gè)時(shí)段均保持著較高的相關(guān)性，決定系數(shù)R2在1967—1979年、1980—2002年和2003—2018年分別為0.96、0.92和0.84。計(jì)算值與實(shí)測(cè)值在1967—2002年和2003—2018年相對(duì)誤差的平均值分別為0.83%和0.45%，計(jì)算效果較好。由圖6(b)可知計(jì)算值與實(shí)測(cè)值在1967—2002年和2003—2018年均基本保持一致。

圖6 累計(jì)淤積量計(jì)算值與實(shí)測(cè)值對(duì)比

伴隨著水庫(kù)控制方式的調(diào)整和來(lái)水來(lái)沙條件的變化，三門峽庫(kù)區(qū)共經(jīng)歷了1960—1969年的快速淤積期、1970—1973年的快速?zèng)_刷期、1974—2002年的緩慢淤積期與2003年之后的緩慢沖刷期4個(gè)階段[23]。4個(gè)時(shí)段中，2003—2018年淤積量的變化范圍最小，為29.41億～30.66億m3，最大、最小值差值僅為1.25億m3。從圖6(b)中可以看出2003—2018年數(shù)據(jù)點(diǎn)(紅點(diǎn))的分布范圍明顯小于1967—2002年數(shù)據(jù)點(diǎn)(藍(lán)點(diǎn))。如何準(zhǔn)確模擬2003年之后累計(jì)淤積量的小幅度變化是滯后響應(yīng)模型研究中的難點(diǎn)，本文通過(guò)改進(jìn)調(diào)整速率參數(shù)的計(jì)算公式提升了模型在2003年之后累計(jì)淤積量的計(jì)算精度，將決定系數(shù)R2從之前的最好計(jì)算結(jié)果0.77[17]升高至目前的0.84。

3.2 平衡值公式起始高程和比降變化范圍及合理性檢驗(yàn)在收集1967—2018年黃淤1—41斷面實(shí)測(cè)汛后深泓高程數(shù)據(jù)的基礎(chǔ)上，得到了黃淤1—41斷面深泓高程的線性擬合公式，并對(duì)汛后深泓縱剖面比降、深泓起始高程以及汛后潼關(guān)站的深泓高程進(jìn)行了統(tǒng)計(jì)，結(jié)果如表1所示。其中，1975—2018年為潼關(guān)高程的計(jì)算時(shí)段，1967—2018年為累計(jì)淤積量的計(jì)算時(shí)段。

1975—2018年汛后潼關(guān)高程的多年平均值為327.55 m，汛后潼關(guān)深泓高程的多年平均值為324.62 m，兩者之差Δh為2.93 m。1975—2018年實(shí)測(cè)深泓線的起始高程變化范圍為291.78～298.79 m，加上Δh，即為式(10)所示Zw的合理變化范圍：

294.71 m≤Zw≤301.72 m

(27)

表1 三門峽水庫(kù)實(shí)測(cè)汛后縱剖面參數(shù)統(tǒng)計(jì)結(jié)果

模型計(jì)算結(jié)果顯示，1975—2018年期間，水庫(kù)沖淤平衡時(shí)深泓縱剖面比降Jes的變化范圍為0.250‰～0.262‰，在同時(shí)段汛后深泓比降Js的變化范圍之內(nèi)(表1)。Jes的時(shí)段平均值為0.255‰，與同一時(shí)段的汛后深泓比降Js多年平均值為0.254‰接近。模型計(jì)算得到Zw的變化范圍為296.80～300.26 m，在式(27)所示的合理變化范圍之內(nèi)。

在累計(jì)淤積量的模擬中，由于概化河底高程與實(shí)測(cè)深泓高程存在差別，所以僅對(duì)概化河底縱剖面平衡比降的變化范圍進(jìn)行討論。

(28)

1967—2018年實(shí)測(cè)深泓線比降的變化范圍為0.144‰～0.326‰，多年平均值為0.253‰。考慮到在平衡條件下比降的變化范圍較小，且為了保證概化河寬B率定結(jié)果的準(zhǔn)確性。使Jem在實(shí)測(cè)深泓線比降多年平均值0.253‰上下約0.022‰的范圍內(nèi)變化：

0.231‰≤Jem≤0.275‰

(29)

根據(jù)實(shí)測(cè)斷面統(tǒng)計(jì)，在1967—2016年，黃淤1—41斷面沖淤河寬時(shí)段均值的變化范圍為0.4～2.3 km，據(jù)此可對(duì)概化河寬B的變化范圍進(jìn)行如下約束：

0.4 km≤B≤2.3 km

(30)

模型率定得到概化河寬B=1.6 km，模型計(jì)算得到Jem的變化范圍為0.237‰～0.273‰，在式(29)所規(guī)定的合理變化范圍之內(nèi)；時(shí)段平均值為0.243‰，同一時(shí)段汛后深泓比降Js的多年平均值為0.253‰，可以發(fā)現(xiàn)Jem和Js的時(shí)段均值接近。因此，本模型中比降Jem的變化范圍是合理的。

4 關(guān)于滯后調(diào)整過(guò)程的討論

4.1 達(dá)到準(zhǔn)平衡所需時(shí)間的計(jì)算在滯后響應(yīng)模型概化示意圖1中，根據(jù)“擾動(dòng)”和“準(zhǔn)平衡態(tài)”兩個(gè)臨界點(diǎn)，將特征變量變化過(guò)程分為擾動(dòng)前、調(diào)整階段、準(zhǔn)平衡階段三個(gè)時(shí)段。據(jù)式(1)可將調(diào)整時(shí)間t表示為：

(31)

將河道受到擾動(dòng)后特征變量從初始值調(diào)整到平衡值所需的變化量Δy稱為目標(biāo)變化量，Δy=(ye-y0)。因?yàn)槠胶庵祔e只可以無(wú)限逼近而無(wú)法達(dá)到，所以還需要定義一個(gè)準(zhǔn)平衡態(tài)。設(shè)定特征變量的變化量達(dá)到目標(biāo)變化量Δy的95%時(shí)河流系統(tǒng)進(jìn)入準(zhǔn)平衡態(tài)，把特征變量y從初始值y0達(dá)到準(zhǔn)平衡態(tài)y′e所需的時(shí)間T稱為準(zhǔn)平衡時(shí)間。即當(dāng)y=y′e時(shí)，y′e-y0=0.95(ye-y0)或y′e=0.95ye+0.05y0，將y=y′e代入式(31)，可得：

(32)

由式(32)計(jì)算得到潼關(guān)高程和累計(jì)淤積量在不同年份的準(zhǔn)平衡時(shí)間如圖7所示。注意到潼關(guān)高程與累計(jì)淤積量準(zhǔn)平衡時(shí)間均在2014年和2015年出現(xiàn)峰值，原因是這兩年的入庫(kù)含沙量嚴(yán)重偏低，分別只有3.18和2.70 kg/m3，分別按式(22)和(25)中βtg淤和βV淤的公式，計(jì)算所得2014年和2015年的潼關(guān)高程調(diào)整速率參數(shù)βtg分別只有0.22 a-1和0.20 a-1(圖8(a))，累計(jì)淤積量調(diào)整速率參數(shù)βv分別只有0.09 a-1和0.06 a-1(圖8(b))。根據(jù)式(32)可知，特征變量的準(zhǔn)平衡時(shí)間與調(diào)整速率參數(shù)β成反比，所以潼關(guān)高程和累計(jì)淤積量準(zhǔn)平衡時(shí)間均在2014年和2015年出現(xiàn)高峰。由此可見(jiàn)，嚴(yán)重偏小的含沙量是導(dǎo)致2014年和2015年調(diào)整時(shí)間偏長(zhǎng)的原因。

圖7 潼關(guān)高程和累計(jì)淤積量準(zhǔn)平衡時(shí)間逐年分布

圖8 潼關(guān)高程和累計(jì)淤積量調(diào)整速率參數(shù)逐年變化

由于2014年和2015年潼關(guān)高程和累計(jì)淤積量準(zhǔn)平衡時(shí)間與前后年份差異較大，為避免個(gè)別年份的極端值掩蓋準(zhǔn)平衡時(shí)間的整體特性，在準(zhǔn)平衡時(shí)間的統(tǒng)計(jì)中略去2014年和2015年的數(shù)據(jù)。同時(shí)，為了便于對(duì)比，將潼關(guān)高程與累計(jì)淤積量準(zhǔn)平衡時(shí)間的研究時(shí)段取為一致，統(tǒng)計(jì)結(jié)果見(jiàn)表2。

表2 潼關(guān)高程與累計(jì)淤積量準(zhǔn)平衡時(shí)間 (單位：a)

統(tǒng)計(jì)結(jié)果顯示，潼關(guān)高程準(zhǔn)平衡時(shí)間多年平均值約為10.2 a，累計(jì)淤積量準(zhǔn)平衡時(shí)間多年平均值約為6.7 a，潼關(guān)高程的準(zhǔn)平衡時(shí)間長(zhǎng)于累計(jì)淤積量的準(zhǔn)平衡時(shí)間表明潼關(guān)高程對(duì)水沙條件的響應(yīng)要慢于累計(jì)淤積量。同時(shí)，注意到2003年之后，受水庫(kù)調(diào)控方式改變及來(lái)沙量減小等因素的影響，潼關(guān)高程和累計(jì)淤積量的準(zhǔn)平衡時(shí)間均有一定程度的延長(zhǎng)。

如果水沙條件不變，自河道受擾動(dòng)開(kāi)始，經(jīng)過(guò)時(shí)間T(準(zhǔn)平衡時(shí)間)后，河道將進(jìn)入準(zhǔn)平衡態(tài)。但由于水沙條件總是處于不斷變化之中，這一準(zhǔn)平衡態(tài)一般難以達(dá)到。此外，庫(kù)區(qū)河道達(dá)到淤積平衡與水庫(kù)淤滿并不是一個(gè)概念，淤積平衡是指一定水沙條件所對(duì)應(yīng)的庫(kù)區(qū)河道平衡形態(tài)，而水庫(kù)淤滿或達(dá)到淤積使用壽命則是指淤沙庫(kù)容淤滿而沒(méi)有足夠的庫(kù)容來(lái)攔截泥沙。在三門峽水庫(kù)蓄清排渾運(yùn)用下，若水沙條件不變，達(dá)到淤積平衡同時(shí)保持一定的淤沙庫(kù)容是可以實(shí)現(xiàn)的。

4.2 特征變量逐時(shí)段以平衡值為目標(biāo)的調(diào)整過(guò)程由于外部條件的變化，特征變量每個(gè)時(shí)段的平衡值yei會(huì)隨著外部條件的改變而處于動(dòng)態(tài)調(diào)整之中。

根據(jù)滯后響應(yīng)模型單步遞推公式(式(17))，可將每個(gè)計(jì)算時(shí)段內(nèi)特征變量平衡值與計(jì)算值的關(guān)系表示為：

(33)

式中特征變量計(jì)算值與平衡值的關(guān)系可進(jìn)一步表達(dá)為：

yi-1

(34)

yi-1>yi>yei(沖刷時(shí))

(35)

圖9(a)(b)分別為潼關(guān)高程和累計(jì)淤積量平衡值與計(jì)算值的關(guān)系，可以看到在模型計(jì)算中，特征變量平衡值起到了“調(diào)整目標(biāo)”的作用。

圖9 特征變量平衡值與計(jì)算值關(guān)系

在每個(gè)計(jì)算時(shí)段內(nèi)，特征變量會(huì)以時(shí)段初值yi-1(即上一時(shí)段末特征變量的計(jì)算值)為基礎(chǔ)，以本時(shí)段的平衡值yei為目標(biāo)，隨著時(shí)間推進(jìn)向yei靠近，經(jīng)過(guò)一個(gè)時(shí)段的調(diào)整得到時(shí)段末計(jì)算值yi。當(dāng)yei高于yi-1時(shí)，特征變量向增大方向調(diào)整；反之，特征變量向減小方向調(diào)整。按照這一模式，特征變量的計(jì)算值從初始年份開(kāi)始，通過(guò)逐時(shí)段調(diào)整得到當(dāng)前年的計(jì)算值。這也說(shuō)明了當(dāng)前時(shí)段的河床演變不僅受當(dāng)前水沙條件的影響，而且通過(guò)邊界條件，還受前期水沙條件的影響[2]。

5 結(jié)論

(1)基于水庫(kù)沖淤平衡縱剖面幾何關(guān)系和機(jī)理分析，推導(dǎo)了潼關(guān)高程和累計(jì)淤積量平衡值的計(jì)算公式；考慮調(diào)整速率隨庫(kù)區(qū)沖淤和水沙條件的變化，提出了滯后響應(yīng)模型調(diào)整參數(shù)β的計(jì)算方法；進(jìn)而采用滯后響應(yīng)模型的單步解析模式，建立了潼關(guān)高程和累計(jì)淤積量的滯后響應(yīng)計(jì)算方法，論證了模型的合理性，計(jì)算結(jié)果表明，潼關(guān)高程1975—2002年和2003—2018年計(jì)算值與實(shí)測(cè)值決定系數(shù)R2分別達(dá)到0.80和0.75，累計(jì)淤積量1967—2002年和2003—2018年計(jì)算值與實(shí)測(cè)值決定系數(shù)R2分別達(dá)到0.95和0.84。(2)定義特征變量的變化量達(dá)到目標(biāo)變化量的95%時(shí)河流系統(tǒng)進(jìn)入準(zhǔn)平衡態(tài)，研究了潼關(guān)高程和累計(jì)淤積量達(dá)到準(zhǔn)平衡所需時(shí)間的變化規(guī)律，結(jié)果顯示，潼關(guān)高程達(dá)到準(zhǔn)平衡的平均時(shí)間約為10.2年，累計(jì)淤積量準(zhǔn)平衡時(shí)間約為6.7年。潼關(guān)高程的準(zhǔn)平衡時(shí)間長(zhǎng)于累計(jì)淤積量的準(zhǔn)平衡時(shí)間，反映了潼關(guān)高程的變化滯后于累計(jì)淤積量的變化。(3)探討了潼關(guān)高程和累計(jì)淤積量逐時(shí)段以平衡值為目標(biāo)的調(diào)整過(guò)程與特點(diǎn)，在每個(gè)時(shí)段內(nèi)，特征變量會(huì)以時(shí)段初值為基礎(chǔ)，以本時(shí)段平衡值為調(diào)整目標(biāo)，隨著時(shí)間推進(jìn)向時(shí)段平衡值靠近。一個(gè)時(shí)段的調(diào)整結(jié)果會(huì)作為下一個(gè)時(shí)段的初始條件使前期的水沙條件對(duì)后期的河床沖淤產(chǎn)生影響，顯示了河床演變的自動(dòng)調(diào)整過(guò)程與滯后響應(yīng)特性。