基于Budyko假設(shè)的黃河下游非一致性徑流模擬

叢培月， 牟獻(xiàn)友， 冀鴻蘭， 張寶森

(1.內(nèi)蒙古農(nóng)業(yè)大學(xué) 水利與土木建筑工程學(xué)院， 內(nèi)蒙古 呼和浩特 010018； 2.黃河水利委員會 黃河水利科學(xué)研究院， 河南 鄭州 450003)

1 研究背景

黃河是我國第二大河流，多年平均徑流量為580×108m3，約占全國總徑流量的2%[1]。 2014年，IPCC全球氣候變化第五次評估報(bào)告顯示：全球氣候變暖已被廣泛認(rèn)可，變暖體現(xiàn)在地表溫度上升、極端氣候事件頻率增加等方面，1880-2012年全球平均溫度已升高0.85℃，氣候變化的嚴(yán)重性遠(yuǎn)超于以往的認(rèn)知[2]。氣候變化會改變降水及蒸散發(fā)，造成水文氣象序列出現(xiàn)非一致性[3]。水文氣象的非一致性會影響陸地水循環(huán)過程，驅(qū)動徑流量的變化[4]。由于歷史資料的缺失，黃河徑流量演化規(guī)律的分析和判別將受到一定的限制[5]，因此，流域徑流關(guān)系演變分析、徑流模擬應(yīng)運(yùn)而生，成為當(dāng)前水文水資源領(lǐng)域研究的前沿與熱點(diǎn)之一，也是國際上一個具有挑戰(zhàn)性的問題[6]。

近50年來，黃河徑流量較以往有著顯著差異，國內(nèi)外許多學(xué)者針對中國各江河流域的徑流模擬，進(jìn)行了定量分析研究[7-8]。顧朝軍等[9]分析了黃河陜縣站1919-1997年天然徑流量變化，指出黃河天然徑流量年際變化存在明顯的階段性及周期性特征；王雁等[10]通過對徑流和氣象觀測數(shù)據(jù)的統(tǒng)計(jì)分析，發(fā)現(xiàn)黃河流域由于氣候變化引起了徑流量的減少。張建云等[11]提出，自20世紀(jì)70年代以來，受氣候變化的影響，黃河中游河川徑流呈現(xiàn)明顯的下降趨勢。孫福寶等[12]基于Budyko假設(shè)的傅報(bào)璞公式檢驗(yàn)了干旱區(qū)水熱耦合理論。熊立華[13]分析了Budyko假設(shè)在中國南部濕潤地區(qū)的適用性，并提出Budyko公式改進(jìn)方法，提高了徑流模擬精度。張麗梅等[14]利用Budyko假設(shè)分析渭河年徑流量變化特征。對徑流的模擬大多是基于水文氣象序列概率分布及其參數(shù)隨時間不變的假定[15]，張波等[16]指出，在環(huán)境變化的大背景下，一致性假設(shè)的水文理論無法精準(zhǔn)地揭示水資源演變規(guī)律，在此背景下提出了非一致性年徑流模擬方法，使結(jié)果更加真實(shí)、可靠；唐雄朋等[17]探討了雅魯藏布江拉孜以上流域非一致性徑流過程模擬與預(yù)測，徑流模擬精度較高；鄧曉宇等[18]對氣候變化和人類活動綜合影響下的撫河流域進(jìn)行徑流模擬研究，發(fā)現(xiàn)納西系數(shù)達(dá)到0.82以上。Milly等[19]指出，在氣候變化和人類活動的影響下，年徑流序列通常是非一致性的，導(dǎo)致許多現(xiàn)有方法徑流模擬的精度較低；Merz等[20]根據(jù)奧地利273個盆地41 a的歷史數(shù)據(jù)，通過校準(zhǔn)模型參數(shù)來進(jìn)行非一致性徑流模擬，提高了徑流模擬精度；Pathiraja等[21]通過使用集合卡爾曼濾波器同時估計(jì)模型參數(shù)，對西澳大利亞實(shí)驗(yàn)區(qū)進(jìn)行徑流預(yù)測，得出預(yù)測結(jié)果可被接受。

目前，對基于Budyko公式改進(jìn)的Zhang模型在黃河下游流域的適用性研究尚有空缺。本文利用灰色關(guān)聯(lián)法[22]與滑動樣本熵等方法對黃河下游水文氣象序列進(jìn)行非一致性識別，將改進(jìn)的Zhang模型應(yīng)用于黃河下游，研究其在非一致性水文氣象序列中模擬徑流的適用性。

2 資料與方法

2.1 研究區(qū)概況

河南鄭州桃花峪以下的黃河河段為黃河下游，地處中國東部，流經(jīng)河南、山東，于山東墾利縣注入渤海。黃河下游段位于東經(jīng)110°～120°、北緯30°～40°之間，河長約786 km，為典型游蕩型河段，河道寬淺，流域面積約2.3×104km2，占全流域面積的3%左右，總落差93.6 m，平均比降0.12%。黃河下游屬于溫帶季風(fēng)氣候，夏季高溫多雨(夏季季平均降水量為371.76 mm，夏季季平均氣溫為26.39℃)，冬季寒冷干燥(冬季季平均降水量為56.71 mm，冬季季平均氣溫為1.34℃)，多年平均降雨量約617.52 mm，年內(nèi)降水分布不均，多年平均氣溫約14.49℃，多年平均徑流量約253.45×108m3，多年平均日照時數(shù)約2 548.77 h。根據(jù)降雨徑流分帶，大部分為半濕潤地帶。

2.2 數(shù)據(jù)來源

收集了黃河下游 7 個水文站(花園口、夾河灘、高村、孫口、艾山、濼口、利津)1979-2015年共37 a的逐日降雨、逐日氣溫、逐日流量等水文資料以及周邊4個氣象站(鄭州、開封、濟(jì)南、惠民)的1979-2015年逐日日照時數(shù)等氣象資料，數(shù)據(jù)均來源于黃河流域歷年水文年鑒、中國氣象數(shù)據(jù)網(wǎng)(http://data.cma.cn)。

2.3 研究方法

本文利用灰色關(guān)聯(lián)法分析徑流復(fù)雜現(xiàn)象與降水、氣溫發(fā)展態(tài)勢的關(guān)聯(lián)、接近的程度，從而確定徑流形成的主控因素。由于水文氣象序列的非一致性會直接影響降水和潛在蒸散量，造成氣候-徑流關(guān)系的方式不同于過去的觀察[23]，因此在非一致性水文氣象序列下，選用基于歷史經(jīng)驗(yàn)的水文模型將影響模擬精度。本文首先進(jìn)行水文氣象序列的一致性驗(yàn)證，利用累積距平法和滑動平均法進(jìn)行初步檢驗(yàn)，并基于滑動樣本熵法[24]進(jìn)行突變分析，聯(lián)合多種分析方法對水文氣象一致性進(jìn)行診斷。由于實(shí)際蒸發(fā)量難以通過觀測直接獲得，通常利用潛在蒸散發(fā)量來計(jì)算實(shí)際蒸發(fā)量[25]，本文采用在濕潤與半濕潤區(qū)優(yōu)先選擇的Hargreaves[26]法估算潛在蒸發(fā)量，通過潛在蒸發(fā)量利用Budyko假設(shè)計(jì)算蒸發(fā)量，最后利用水量平衡法反算徑流量。

2.3.1 Budyko假設(shè) Budyko(1974年)假設(shè)基本原理是：假定流域內(nèi)降水和潛在蒸散發(fā)量決定了多年平均蒸散發(fā)，蒸散發(fā)量為降水與潛在蒸散發(fā)量之間的函數(shù)，即可以對降水、潛在蒸散發(fā)和徑流之間定量關(guān)系進(jìn)行較好的擬合[27]。但是由于植被類別、土壤水力特性及地形的區(qū)別，使得Budyko假設(shè)的廣泛性受到限制，Zhang等在Budyko基礎(chǔ)上引進(jìn)模型參數(shù)ω，改進(jìn)了Budyko假設(shè)公式。

Budyko假設(shè)可簡單表達(dá)為：

(1)

(2)

Zhang公式：

(3)

式中：E為年平均蒸發(fā)量，mm；P為年平均降水量，mm；ω為植被水分利用系數(shù)；F為干燥指數(shù)；E0為年平均潛在蒸散發(fā)量，mm。

水量平衡法：

R=P-E

(4)

式中：R為徑流深，mm；P為降水量，mm；E為蒸發(fā)量，mm。

2.3.2 模擬精度分析 本文采用納西系數(shù)(NSE)對年徑流量的模擬精度進(jìn)行評價。NSE∈(-∞,1)，其中當(dāng)NSE=1時，模擬精度最優(yōu)；NSE∈(0.01,1)時，模擬精度為可接受。

(5)

3 結(jié)果與分析

3.1 降水、氣溫與徑流灰色關(guān)聯(lián)度

以黃河下游7個水文站為代表站，設(shè)定降雨徑流灰色關(guān)聯(lián)度為γ1，氣溫徑流灰色關(guān)聯(lián)度為γ2。通過計(jì)算得出1979-2015年7個水文站的γ1、γ2，如表1所示。7個站的γ1、γ2均大于0.5，表示徑流量與降水和氣溫關(guān)聯(lián)度較高，即氣溫和降水均為黃河下游徑流量的主控因素。此結(jié)果與張國宏等[28]的研究結(jié)論相似，黃河流域徑流量與氣溫和降水均有較好的相關(guān)性，其中各站均為γ2>γ1，表明在黃河下游徑流形成過程中，氣溫與徑流量的關(guān)聯(lián)度較大，對徑流的影響所占的權(quán)重大于降水。

表1 1979-2015年降水、氣溫與徑流的灰色關(guān)聯(lián)度

3.2 降水、氣溫變化趨勢分析

表2為1979-2015年黃河下游7個水文站的降水與氣溫突變年份。由表2可知，1979-2015年黃河下游流域降水和氣溫突變點(diǎn)不一致，如夾河灘站降水時間序列發(fā)生突變的時間為1999年，氣溫時間序列發(fā)生突變的時間為1993年；突變時間均在20世紀(jì)90年代前后，水文氣象時間序列突變時間具有較高的同步性。本文以花園口站降水、氣溫變化趨勢分析為例，其余6個站不再贅述。

表2 1979-2015年黃河下游7個站降水與氣溫突變年份

黃河下游1979-2015年間多年平均降水量為617.52 mm，降水量總體偏少。由線性擬合結(jié)果(圖1(a))，降水量傾向率為-4.5 mm/10a，呈顯著性減少趨勢(p<0.05)，其中降水量最大值出現(xiàn)在1983年(962.94 mm)，最小值出現(xiàn)在1997年(368.96 mm)，極值比1.61。從降水量累積距平曲線(圖1(b))可以看出：花園口站在1979-2015年間降水量變化呈現(xiàn)出先波動下降、后波動上升的趨勢，即降水量變化分為兩個階段：第1階段1979-2002年，該時期降水量平均值為596.24 mm；第2階段2003-2015年，該時期降水量平均值為639.39 mm，第1階段到第2階段的變化呈現(xiàn)了降水枯-豐轉(zhuǎn)變的過程，表明降水時間序列已發(fā)生變化。

1979-2015年，花園口站多年平均氣溫為15.14℃，年平均氣溫傾向率為0.6℃/10a。由線性擬合(圖2(a))可看出，花園口1979-2015年38 a的年平均氣溫呈現(xiàn)顯著上升趨勢(p<0.05)。其中最高值和最低值分別出現(xiàn)在2015和1984年，分別為16.67和13.68℃，極值比為0.22。與降水量的變化相比，氣溫的變化較為平緩。從圖2(b)的累積距平曲線可以看出，花園口站年平均氣溫表現(xiàn)出先下降后上升的發(fā)展趨勢，自1997年起，花園口站氣溫開始由低溫期向高溫期過渡，表明氣溫時間序列也已發(fā)生變化。

圖1 花園口站1979-2015年年降水量時間序列變化趨勢

圖2 花園口站1979-2015年年均氣溫時間序列變化趨勢

降水量的滑動樣本熵(SampEn)呈現(xiàn)逐漸增加的趨勢(圖3(a))，降水量變化分為3個階段： 第1階段(1979-1986年)、第2階段(1987-1993年)、第3階段(1994-2015年)，即發(fā)生了兩次突變，對應(yīng)的突變點(diǎn)年份分別為1986和1993年，突變點(diǎn)前后的SampEn序列差異性顯著。氣溫的SampEn值在0.6～0.9之間(圖3(b)，其變化也呈現(xiàn)出3個階段：第1階段(1979-1991年)、第2階段(1992-1994年)、第3階段(1995-2015年)，對應(yīng)的突變點(diǎn)年份分別為1991和1994年，氣溫時間序列的SampEn值均大于降水量序列，表明氣溫系統(tǒng)的復(fù)雜程度高于降水系統(tǒng)。

圖3 花園口站1979-2015年年降水量、氣溫時間序列SampEn分析

整體而言，從發(fā)生突變的年份來看，年降水量和年平均氣溫發(fā)生突變的時間有較高的對應(yīng)關(guān)系。由以上水文序列發(fā)生突變可知，水文序列的統(tǒng)計(jì)規(guī)律并非一致，因此水文序列不滿足一致性。非一致性在徑流、降水和氣溫中已被廣泛認(rèn)定。

3.3 潛在蒸散發(fā)量估算

利用Hargreaves方法在黃河下游4個水文站點(diǎn)(花園口、夾河灘、艾山、利津)中進(jìn)行潛在蒸散發(fā)量估算，估算結(jié)果見圖4。從圖4可知，具有相似的規(guī)律特征，年平均潛在蒸散發(fā)量從850～1 200 mm不等，潛在蒸散發(fā)量均呈現(xiàn)隨時間減少的顯著趨勢(p<0.05)，4個站點(diǎn)的氣候傾向率分別為-20.8、-34.8、-76.3和-38.2 mm/10a，其中利津站在1990年前后出現(xiàn)了小幅回升，1979-1993年潛在蒸散發(fā)量呈現(xiàn)下降的趨勢，氣候傾向率為-34.60 mm/10a,2010-2014年潛在蒸散發(fā)量呈現(xiàn)上升的趨勢，氣候傾向率為51.79 mm/10a，這些變化符合近50年來中國潛在蒸散量的變化趨勢[29]。

3.4 Budyko假設(shè)徑流模擬精度分析

由于Budyko假設(shè)沒有考慮人類活動和氣候變化對徑流量的影響，因此利用降水量與潛在蒸散發(fā)量，同時基于Budyko假設(shè)改進(jìn)的Zhang模型求得蒸發(fā)量，通過水量平衡法反算徑流量。本文以1982年為例，模擬計(jì)算了逐月徑流量并與實(shí)測值進(jìn)行對比，如圖5所示。分析圖5可知：花園口站在8月模擬結(jié)果誤差最大(NSE<0)，在2月(NSE=0.81)和5月(NSE=0.97)模擬結(jié)果較為理想，夾河灘站和艾山站分別在7、8月模擬結(jié)果誤差最大(NSE=0)，利津站在4月(NSE=0.97)、5月(NSE=0.99)、6月(NSE=0.86)模擬結(jié)果較為理想，模擬精度較高。表3給出了黃河下游4個典型區(qū)1982年月均徑流模擬NSE檢驗(yàn)結(jié)果，4個典型區(qū)的NSE均小于0。綜合以上研究結(jié)果，在月尺度徑流模擬中各個月模擬值與實(shí)測值之間誤差較大，模擬精度過低，因此Zhang公式不適用于月尺度徑流模擬中。

利用Zhang公式對4個水文站的1979-2015年徑流過程進(jìn)行模擬，結(jié)果如圖6所示，模擬的徑流變化過程與實(shí)際的徑流變化過程基本吻合，能較好地反映徑流的變化趨勢和數(shù)值，但在2000年后個別年份模擬值與實(shí)際值差別較大，如花園口站2005年(NSE=-0.62)、2011年(NSE=-2.27)，艾山站2010年(NSE=-4.77)。但整體來看，4個典型區(qū)徑流模擬的納西系數(shù)均大于0.9(如表4所示)，說明模型精度較高，模擬結(jié)果可以接受，Zhang公式在年徑流量模擬方面具有一定的適用性。

圖4 黃河下游4個水文站1979-2015年潛在蒸散發(fā)量變化

圖5 1982年黃河下游4個水文站月徑流量觀測值與模擬值變化分析

圖6 1979-2015年黃河下游4個水文站年徑流量觀測值與模擬值變化分析

表3 1982年黃河下游4個水文站月徑流量NSE對比

表4 1979-2015年黃河下游4個水文站年徑流量NSE對比

徑流量模擬的差異性受諸多因素影響，水文氣象序列的選擇、人類活動和氣候變化等可直接影響徑流量模擬精度。本文利用基于Budyko假設(shè)改進(jìn)的Zhang模型對黃河下游進(jìn)行徑流量模擬，經(jīng)過驗(yàn)證表明模型精度較高，適應(yīng)于黃河下游流域非一致性水文氣象序列的徑流量模擬，可為黃河下游徑流量分析、預(yù)報(bào)提供可靠的數(shù)據(jù)支持。

4 結(jié) 論

(1)黃河下游降水量呈現(xiàn)出顯著減少趨勢(p<0.05)，氣溫呈現(xiàn)出顯著增加趨勢(p<0.05)，利用滑動樣本熵從動力學(xué)的角度判別出降水量、氣溫均在20世紀(jì)90年代前后出現(xiàn)突變，且突變年份具有較高的同步性。

(2)黃河下游徑流變化過程受多因素影響，其中γ1和γ2均大于0.5，表明降水量和氣溫與徑流量的關(guān)聯(lián)程度均較高，是影響徑流量變化的主要因素。

(3)基于Budyko假設(shè)改進(jìn)的Zhang模型模擬了黃河下游非一致性水文氣象序列下的年徑流量，NSE達(dá)到0.9以上，模擬效果良好，表明Zhang模型在黃河下游非一致性水文氣象序列中模擬徑流量具有可適性。