張 茹, 張志全, 羅 凱

(1.山西水利職業(yè)技術(shù)學(xué)院, 山西 運(yùn)城 044000; 2.河海大學(xué) 水利水電學(xué)院,南京 210098; 3.黃河水利委員會(huì) 上游水文水資源局, 蘭州 730030)

黃河下游河道的輸沙功能是決定黃河入海泥沙通量的重要因素，對(duì)黃河口海岸的演變、黃海與渤海水下三角洲的沉積具有深遠(yuǎn)影響[1]。在不同時(shí)期，黃河下游河道輸沙功能的演變趨勢(shì)與其造床造陸的作用在自然和人為因素的共同驅(qū)動(dòng)下動(dòng)態(tài)變化[2]。自20世紀(jì)80年代中期以來(lái)，黃河流域降水周期性減少，同時(shí)農(nóng)業(yè)灌溉需求的增大使得河道引水量也大幅增長(zhǎng)，黃河下游河道輸沙功能迅速減弱，河道發(fā)生萎縮，甚至出現(xiàn)季節(jié)性斷流，三角洲造陸功能也顯著退化，引起各方廣泛關(guān)注[3-4]。隨著小浪底水庫(kù)的建成，水庫(kù)的運(yùn)行調(diào)度能夠?qū)S河下游水沙過(guò)程起到一定的調(diào)節(jié)作用，進(jìn)而促使河道輸沙功能發(fā)生與之前截然不同的變化。許炯心[5]就指出，通過(guò)小浪底水庫(kù)進(jìn)行調(diào)水調(diào)沙，可以使得場(chǎng)次洪水平均含沙量維持在35 kg/m3左右的輸沙最優(yōu)含沙量。這種調(diào)水調(diào)沙使得黃河下游河道輸沙功能大幅增強(qiáng)，扭轉(zhuǎn)了黃河下游河道逐年萎縮的態(tài)勢(shì)[6]。近年來(lái)經(jīng)濟(jì)社會(huì)的高速發(fā)展促使人們?cè)邳S河流域的治理與開(kāi)發(fā)活動(dòng)增加，特別是水土流失綜合治理措施的大范圍實(shí)施和黃河下游農(nóng)業(yè)灌溉用水的持續(xù)增大，加之全球氣候變暖，黃河下游河道的水沙條件將會(huì)進(jìn)一步改變，也會(huì)導(dǎo)致其輸沙功能的響應(yīng)與變化[7]。史紅玲等[8]采用M-K檢驗(yàn)法及秩和檢驗(yàn)法對(duì)黃河干流主要水文站水沙序列進(jìn)行分析，發(fā)現(xiàn)水土保持措施減沙及水庫(kù)建設(shè)攔沙對(duì)近年來(lái)黃河來(lái)沙量的持續(xù)減少起到了主導(dǎo)作用。胡春宏等[9]的研究表明，黃河河道輸沙能力與水利水保工程建設(shè)呈現(xiàn)協(xié)同的階段變化，近年來(lái)黃河干流利于輸沙的流量持續(xù)時(shí)間和水量、沙量都大幅減少，需對(duì)水沙變化機(jī)理深入研究，以準(zhǔn)確預(yù)測(cè)新情勢(shì)下黃河水沙變化趨勢(shì)。夏軍等[10]則指出，在氣候變化條件下，應(yīng)完善小浪底水庫(kù)水沙調(diào)控體系，探討黃河下游河道更加高效的輸沙模式。

目前，對(duì)于在強(qiáng)人類(lèi)活動(dòng)影響下黃河下游輸沙量變化及其所引起的河床演變和三角洲造陸變化等問(wèn)題已有諸多研究成果，但對(duì)近年來(lái)黃河下游河道輸沙功能變化及其影響因素的分析研究尚且不足。探明近年來(lái)強(qiáng)人類(lèi)活動(dòng)影響下黃河下游輸沙功能演變的原因及其主導(dǎo)因素，是當(dāng)前黃河下游開(kāi)發(fā)、治理與保護(hù)迫切需要解決的問(wèn)題。本文基于黃河下游河段各站點(diǎn)2007—2020年實(shí)測(cè)徑流與泥沙資料，計(jì)算分析近年黃河下游河道輸沙功能的時(shí)空演變特征，并結(jié)合流域降水、水土保持和引水引沙數(shù)據(jù)，討論各因素對(duì)輸沙功能演變的影響，為黃河下游河道輸沙與河床演變的有關(guān)研究提供新的認(rèn)識(shí)。

1 研究區(qū)概況

黃河發(fā)源于青藏高原巴顏喀拉山脈，流經(jīng)青海、四川等9個(gè)省區(qū)，全長(zhǎng)約5 464 km，其流域面積約為75萬(wàn)km2，是世界上含沙量最大的河流之一。黃河流域主要的泥沙來(lái)源為中游的黃土高原，區(qū)域產(chǎn)沙量約占黃河下游來(lái)沙總量的90%[11]。黃河下游指河南鄭州桃花峪以下的黃河河段，河段平均比降約為0.12‰。本文研究河段為黃河下游花園口—利津河段(小浪底水庫(kù)下游110～735 km)，全長(zhǎng)約625 km。河段中還設(shè)有高村(小浪底水庫(kù)下游283 km處)和艾山(小浪底水庫(kù)下游465 km處)兩個(gè)水沙監(jiān)測(cè)站，將河段分3部分。其中花園口—高村河段全長(zhǎng)約173 km，按河床演變特點(diǎn)區(qū)分屬于游蕩型河段，比降較大，河道寬淺，洲灘眾多，河床由易沖刷的細(xì)沙組成；高村—艾山河段全長(zhǎng)約182 km，屬于過(guò)渡型河段，兩岸土質(zhì)較好，比降介于上游游蕩段與下游彎曲段之間；艾山—利津河段則屬于彎曲型河段，全長(zhǎng)約270 km，比降較小，河岸抗沖性較強(qiáng)，加之人為控導(dǎo)工程的作用，其河勢(shì)較穩(wěn)定。

2 數(shù)據(jù)來(lái)源與研究方法

2.1 水文數(shù)據(jù)

本文所用的花園口(113.7°E,34.9°N)、高村(115.1°E,35.4°N)、艾山(116.3°E,36.3°N)和利津(118.3°E,37.5°N)水文站的徑流輸沙資料、泥沙粒徑和河段引水引沙數(shù)據(jù)來(lái)源于《中國(guó)河流泥沙公報(bào)》(2007—2020年)；區(qū)域(山西省、陜西省、甘肅省、內(nèi)蒙古自治區(qū)和寧夏回族自治區(qū))水土流失綜合治理面積數(shù)據(jù)引自《中國(guó)水土保持公報(bào)》(2007—2020年)；流域降水?dāng)?shù)據(jù)源自《中國(guó)水資源公報(bào)》和《水情年報(bào)》(2007—2020年)。

2.2 河段輸沙功能的計(jì)算方法

河道輸沙功能指某一河道在宏觀意義上輸送泥沙的能力，它不要求滿足輸沙平衡條件，因?yàn)闊o(wú)論處于平衡狀態(tài)或非平衡狀態(tài)(淤積和沖刷)下的河道都具有輸沙功能。因此，以黃河下游各河段為單元，從泥沙收支平衡的概念(Sediment budget)出發(fā)[2,12]，黃河下游河道的輸沙功能指標(biāo)Fs(也有學(xué)者稱(chēng)之為河道排沙比)可按下式進(jìn)行計(jì)算：

(1)

式中:u,v分別為河段的進(jìn)口點(diǎn)和出口點(diǎn);Qs,u和Qs,v分別為河段進(jìn)口處和出口處的輸沙量;由于黃河下游灌溉引沙量較大,不可忽略,故分母處需減去河段區(qū)間引沙量Qdiv,uv。

注意到，若取Qs,u為花園口站輸沙量；Qs,v為利津站輸沙量；Qdiv,uv為花園口—高村、高村—艾山、艾山—利津三段引沙量之和，則可計(jì)算得到黃河下游整體的輸沙功能指標(biāo)Fs。

2.3 相對(duì)權(quán)重法

相對(duì)權(quán)重法(relative weight)是由Johnson等[13]于2000年提出的一種估計(jì)線性模型中自變量相對(duì)重要性的方法。本文使用該方法來(lái)量化各影響因素(來(lái)沙系數(shù)、中值粒徑、降水量、水土流失綜合治理面積和引水量)對(duì)黃河下游河道輸沙功能變化的相對(duì)貢獻(xiàn)率。該方法首先以自變量分解出的正交變量作為中間載體，然后分別通過(guò)自變量對(duì)正交變量、因變量對(duì)正交變量建立線性回歸方程，最后自變量的相對(duì)重要性則為兩組回歸系數(shù)平方的積的和。相對(duì)權(quán)重法的主要步驟簡(jiǎn)述如下：

(1) 將自變量(來(lái)沙系數(shù)、中值粒徑、降水量、水土流失綜合治理面積和引水量)和因變量(黃河下游河道輸沙功能)進(jìn)行標(biāo)準(zhǔn)化處理(X,Y)；

(2) 對(duì)標(biāo)準(zhǔn)化后的自變量矩陣X進(jìn)行奇異值分解，獲得相互獨(dú)立的正交變量Z；

(3) 建立正交變量Z與標(biāo)準(zhǔn)化的自變量X和標(biāo)準(zhǔn)化因變量Y之間的線性回歸方程，分別得到回歸系數(shù)λkj(自變量X與正交變量Z之間)和Bk(因變量Y與正交變量Z之間)；

(4) 由于矩陣Z中的正交變量間是不相關(guān)的，因此避免了變量之間的相互作用(共線性)問(wèn)題?？捎忙薻j2表示正交變量矩陣Z所占原始變量矩陣X的貢獻(xiàn)比例。最后，各變量的權(quán)值εj可表示為：

(2)

可以看到，各自變量的權(quán)值之和等于整個(gè)模型的決定系數(shù)R2。最后，各自變量的相對(duì)權(quán)重ψj(即本研究中的相對(duì)貢獻(xiàn)率)可表示為：

(3)

3 黃河下游河道輸沙功能的時(shí)空變化

圖1統(tǒng)計(jì)了黃河下游各河段輸沙功能指標(biāo)在2007—2020年的變化情況?？傮w來(lái)看，花園口—利津河段輸沙功能呈現(xiàn)出明顯的下降趨勢(shì)，至2018—2020年，該段輸沙功能已降至0.97左右，相比于2007年的2.28降幅達(dá)57%。黃河下游花園口—利津河段是黃河向海輸送泥沙的重要輸沙帶，該段輸沙功能的大幅度退化，是近年黃河三角洲造陸速度減緩的重要原因之一。分段來(lái)看，黃河下游花園口—高村段(游蕩段)的輸沙功能最大、高村—艾山段(過(guò)渡段)次之，而艾山—利津段(彎曲段)則最小。3個(gè)河段輸沙功能在2007—2020年均出現(xiàn)不同程度的下降，其中花園口—高村段下降趨勢(shì)較高村—艾山段和艾山—利津段更為顯著，這與顏明等[14]基于長(zhǎng)序列水沙數(shù)據(jù)得出的高村以上河道沖刷強(qiáng)度大于高村以下河段沖刷強(qiáng)度的認(rèn)識(shí)基本相符。

圖1 黃河下游各河段輸沙功能指標(biāo)的變化趨勢(shì)

4 輸沙功能演變的影響因素分析

4.1 來(lái)沙系數(shù)和泥沙粒徑

圖2是2007—2020年黃河下游花園口、高村、艾山、利津4個(gè)水文站來(lái)沙系數(shù)(定義為ξ=S/Q，S為含沙量，Q為流量)和床沙中值粒徑均值的年際變化情況。可以看到，受小浪底水庫(kù)攔沙作用的影響，黃河下游來(lái)沙系數(shù)呈現(xiàn)出減小的趨勢(shì)。已有研究表明，來(lái)沙系數(shù)越大，則河道輸沙功能指標(biāo)越低[2]。因此黃河下游來(lái)沙系數(shù)的減小理應(yīng)使得其輸沙功能指標(biāo)升高，但由前文分析可知，在研究時(shí)期內(nèi)黃河下游河道輸沙功能指標(biāo)整體上表現(xiàn)為下降的趨勢(shì)，這說(shuō)明來(lái)沙系數(shù)的變化并不是該時(shí)期影響輸沙功能演變的主導(dǎo)因素。

圖2 黃河下游來(lái)沙系數(shù)與床沙中值粒徑的變化趨勢(shì)

觀察黃河下游床沙中值粒徑的變化可以發(fā)現(xiàn)，研究時(shí)期內(nèi)中值粒徑呈現(xiàn)增大的趨勢(shì)，這是小浪底水庫(kù)建成后黃河下游河道持續(xù)沖刷而導(dǎo)致的河床粗化現(xiàn)象[15]。由于河道水流挾沙力與泥沙沉速負(fù)相關(guān)，而泥沙沉速與泥沙粒徑正相關(guān)，因此在同等水力條件下河道輸送粗顆粒泥沙的能力要低于輸送細(xì)泥沙的能力。這表明研究時(shí)期內(nèi)床沙中值粒徑的持續(xù)粗化是該時(shí)期黃河下游河道輸沙功能指標(biāo)下降的原因之一。

4.2 降水量和水土保持

圖3是2007—2020年黃河流域降水量及黃河流域中上游主要產(chǎn)沙區(qū)域(山西省、陜西省、甘肅省、內(nèi)蒙古自治區(qū)和寧夏回族自治區(qū))水土流失綜合治理面積的變化情況。如圖3所示，流域降水量表現(xiàn)為周期性波動(dòng)的趨勢(shì)，在研究時(shí)期內(nèi)無(wú)明顯的上升或下降趨勢(shì)，因此降水量的變化不是該時(shí)期河道輸沙功能下降的主要原因。除降水量外，水土保持是影響黃河流域產(chǎn)沙量的另一重要因素[16]。隨著水土保持措施的大范圍實(shí)施，在2007—2020年，黃河流域中游主要產(chǎn)沙區(qū)的水土流失綜合治理面積呈現(xiàn)出明顯的增大趨勢(shì)，2020年相比于2007年水土流失綜合治理面積增幅超過(guò)7.2%。水土流失綜合治理面積的增大使得流域內(nèi)坡面侵蝕強(qiáng)度降低，泥沙匯入量減小，削減了黃河下游的泥沙來(lái)源，是黃河下游河道輸沙功能發(fā)生變化的驅(qū)動(dòng)因素之一。

圖3 流域降水量及黃河下游水土流失

4.3 農(nóng)業(yè)灌溉用水

由于黃河下游灌溉引水量較大，引水的同時(shí)不可避免地引走河道內(nèi)的泥沙，因此河段的區(qū)間引水勢(shì)必會(huì)對(duì)其輸沙功能產(chǎn)生重要影響。本文統(tǒng)計(jì)了2007—2020年黃河下游花園口—利津河段區(qū)間引水量的變化情況(圖4)，在經(jīng)濟(jì)社會(huì)不斷發(fā)展的大背景下，研究時(shí)期內(nèi)黃河下游區(qū)間引水量也顯著增大。在相同的水流含沙量條件下，引水量的增大意味著引沙量也相應(yīng)增加，河道內(nèi)可輸送的沙量減少。因此，研究時(shí)期內(nèi)引水量的增大也是該時(shí)期黃河下游河道輸沙功能指標(biāo)下降的原因之一。

4.4 量化各影響因素對(duì)于輸沙功能演變的貢獻(xiàn)率

一般而言，來(lái)沙系數(shù)、泥沙粒徑、降水量、水土流失綜合治理面積、引水量與輸沙功能指標(biāo)之間的關(guān)系并不是獨(dú)立的，它們之間存在某種程度上的共線性，如水土流失綜合治理會(huì)減少流域產(chǎn)沙，從而影響來(lái)沙系數(shù)。因此本文將輸沙功能指標(biāo)序列作為因變量，上述5種影響因素作為自變量，運(yùn)用相對(duì)權(quán)重法探討上述影響因素在黃河下游河道輸沙功能演變過(guò)程中的相對(duì)貢獻(xiàn)率。表1的計(jì)算結(jié)果表明，在研究時(shí)期內(nèi)，床沙中值粒徑是影響輸沙功能變化最主要因素，其貢獻(xiàn)率達(dá)到37.6%；其次是水土流失綜合治理面積，其對(duì)輸沙功能變化的影響程度約為21.6%；引水量的影響程度在5個(gè)因素中排序第三，為17.5%；來(lái)沙系數(shù)對(duì)輸沙功能變化的影響程度則較小，為15.1%；而降水量的影響最小，僅為8.2%。取影響程度最大的3個(gè)因素——中值粒徑、水土流失綜合治理面積和引水量，將其貢獻(xiàn)率求和可達(dá)到76.7%，表明這3個(gè)因素的變化對(duì)于黃河下游河道輸沙功能演變起到主導(dǎo)作用。

圖4 黃河下游花園口-利津段區(qū)間引水量的變化趨勢(shì)

表1 各影響因素對(duì)黃河下游河道輸沙功能演變的貢獻(xiàn)度

為進(jìn)一步探明黃河下游河道輸沙功能變化與其主要驅(qū)動(dòng)因素之間的線性關(guān)系，構(gòu)建包含泥沙粒徑、水土流失綜合治理面積和引水量這3個(gè)自變量和輸沙功能這一因變量之間的線性回歸模型，基于最小二乘法擬合得到回歸方程表達(dá)式如下：

Fs,uv=-98.016x1+0.104x2-0.01x3+1.196

(4)

式中:x1,x2,x3分別為泥沙粒徑、水土流失綜合治理面積和引水量。由回歸方程的系數(shù)可知,床沙中值粒徑和引水量與河道輸沙功能呈負(fù)相關(guān),與前文的定性分析相符;而水土流失綜合治理面積與河道輸沙功能呈正相關(guān),這是因?yàn)樗帘３执胧┦沟命S河流域中游產(chǎn)沙減小,導(dǎo)致含沙量較小的水流進(jìn)入黃河下游,使得下游河道輸沙能力增強(qiáng)。2007—2020年,黃河下游河道輸沙功能在床沙粗化、水土流失治理面積增大和引水量增加3個(gè)主導(dǎo)因素的影響下,整體上呈現(xiàn)出下降的趨勢(shì)。

為判別回歸方程擬合質(zhì)量的好壞，計(jì)算回歸方程計(jì)算值與實(shí)測(cè)值的均方根誤差(RMSE)和納什效率系數(shù)(NSE)如下：

(5)

(6)

如圖5所示，將回歸方程計(jì)算得到的河道輸沙功能與其實(shí)測(cè)值進(jìn)行對(duì)比，我們發(fā)現(xiàn)計(jì)算值與實(shí)測(cè)值基本吻合，表明考慮泥沙粒徑、水土流失綜合治理面積和引水量這3個(gè)主要影響因素的多元回歸方程可以較為準(zhǔn)確地反映輸沙功能的變化情況。

圖5 輸沙功能回歸方程計(jì)算值與實(shí)測(cè)值的對(duì)比

5 結(jié) 論

(1) 黃河下游河道輸沙功能指標(biāo)在2007—2020年呈現(xiàn)顯著下降的趨勢(shì)，研究時(shí)期內(nèi)降幅達(dá)57%。花園口—高村段輸沙功能比高村—艾山段和艾山—利津段要大，其降幅也更大。黃河下游河道輸沙功能的大幅度退化，是近年黃河三角洲造陸速度減緩的重要原因之一。

(2) 2007—2020年造成黃河下游河道輸沙功能下降的主導(dǎo)因素是河床粗化、水土流失綜合治理面積增加和引水量增大，三者對(duì)輸沙功能變化的相對(duì)貢獻(xiàn)率分別達(dá)到37.6%，21.6%，17.5%。其中，河床粗化與引水量增大使得河道輸沙功能減弱，而水土流失綜合治理面積增加使得河道輸沙功能增強(qiáng)。此外，通過(guò)多元線性回歸，考慮上述3個(gè)主導(dǎo)因素得到的回歸方程可以較為準(zhǔn)確地反映輸沙功能的年際變化情況。