潘佑靜 唐麗霞 楊 智 杜 波 羅 金

(1. 貴州省甕安縣猴場水庫管理所，貴州 黔南 550400；2. 貴州大學(xué)林學(xué)院，貴州 貴陽 550025；3. 貴州省水土保持技術(shù)咨詢研究中心，貴州 貴陽 550002；4. 貴州省清鎮(zhèn)水務(wù)局，貴州 貴陽 551400；5. 貴州省甕安縣天文水庫管理所，貴州 甕安 550400)

喀斯特地區(qū)是典型的生態(tài)環(huán)境脆弱區(qū)，下墊面的高度異質(zhì)性加之人為破壞使得該地區(qū)的水土流失極為嚴(yán)重，土壤侵蝕產(chǎn)沙的成因非常復(fù)雜[1]。長期的巖溶作用構(gòu)成地表和地下雙層空間結(jié)構(gòu)[2]，降雨是造成該地區(qū)地表土壤水蝕的主要動力因素，如何確定有效的侵蝕性降雨及其對地表侵蝕產(chǎn)沙的影響尤為重要[3]。目前對侵蝕性降雨的研究主要集中在雨強(qiáng)、降雨歷時及降雨量等方面[3-4]，對前期有效降雨研究較少考慮。前期有效降雨通過影響坡面土壤的含水量，改變土壤的粘性、抗剪強(qiáng)度、入滲初期的水勢梯度等物理特性，影響土壤水文過程及土壤侵蝕。國內(nèi)外對前期有效降雨的研究主要集中在有效前期降雨值的確定[5-7]，泥石流形成中的作用[8]，以及常態(tài)地貌下與流域產(chǎn)流產(chǎn)沙的關(guān)系研究[9]，相對常態(tài)地貌，前期有效降雨對喀斯特地區(qū)地表產(chǎn)流產(chǎn)沙影響更為顯著，然而這方面現(xiàn)有的研究卻很少?？λ固氐貐^(qū)碳酸鹽巖風(fēng)化發(fā)育的石灰土土層薄、土壤總量少、儲水能力低，流域以蓄滿產(chǎn)流為基本產(chǎn)流機(jī)制[10]，在非超滲產(chǎn)流機(jī)制主導(dǎo)下，前期土壤水分對入滲能力有著直接的影響，降雨發(fā)生時的土壤水分越高，土壤入滲能力就越弱，從而更有利于地表徑流的產(chǎn)生[11]，產(chǎn)沙量也隨之呈現(xiàn)增加的變化趨勢[12-13]。

本研究以貴州省喀斯特地區(qū)的關(guān)嶺縣享樂小流域為研究區(qū)，通過實測的日降雨量，采用產(chǎn)沙累積頻度-降雨量分形關(guān)系計算流域前期有效降雨量[14-15]，劃分流域發(fā)生或不發(fā)生侵蝕產(chǎn)沙的降雨條件，從而獲取喀斯特區(qū)地表侵蝕產(chǎn)沙的降雨條件閾值，以期為貴州省喀斯特地區(qū)防治坡面水土流失、坡耕地整治、生態(tài)恢復(fù)與重建提供參考依據(jù)。

1 研究區(qū)概況

研究區(qū)為貴州省關(guān)嶺縣的享樂小流域，地處東經(jīng)105°32′59″～105°34′34″，北緯25°48′22″～25°50′12″，流域面積為3.93 km2。享樂小流域?qū)匍L江流域烏江水系，地貌類型是以溶峰林、峰叢-洼地、谷地為主的巖溶低中山丘陵，流域內(nèi)有地下水出露點與落水洞，流域巖溶面積占83.84%，水土流失以水力侵蝕為主。地勢為西高東低，由西北向東南傾斜，西部最高達(dá)1 775 m，中部次之為1 470 m，東部最低點1 160 m，相對高差為615 m。流域地處亞熱帶濕潤季風(fēng)氣候區(qū)，年平均氣溫19.2 ℃；月平均最低氣溫1月最低，為6.5 ℃； 8月最高，為24.1 ℃。年平均降雨量為1 236 mm，無霜期298 d。土壤主要以水稻土、黃壤、石灰土為主；現(xiàn)有植被類型主要為次生喬木林、疏幼林，森林覆蓋率為22.42%，樹種以杉木 (Cunninghamialanceolata)、柳杉 (Cryptomeriafortunei)、光皮樺 (Betulaluminifera) 為主。

2 研究方法

2.1 數(shù)據(jù)觀測

于流域小型氣象站和卡口站內(nèi)進(jìn)行觀測，采用HOBO小型自動氣象站記錄流域歷次降雨過程數(shù)據(jù)，包括降雨量，降雨開始及停止時間，降雨歷時，共獲取流域2009—2013年5 a間的降雨事件數(shù)據(jù)?？谡静捎镁匦伪”诹克?，徑流數(shù)據(jù)采用SWZ 型周記式水位計進(jìn)行記錄，同時采用自動水位取樣器取樣，將水樣用濾紙過濾，獲取河道泥沙懸移質(zhì)、推移質(zhì)數(shù)據(jù)，計算出每次降雨事件中流域產(chǎn)沙量，2009—2013年共有562 d發(fā)生了降雨，其中70 d的降雨造成了侵蝕產(chǎn)沙，詳細(xì)數(shù)據(jù)見表1。

表1 2009—2013年降雨、產(chǎn)沙天數(shù)統(tǒng)計Table 1 The statistics of rain fall and sediment during 2009-2013

2.2 分形分布擬合

分形分布廣泛存在于自然界，所有規(guī)模的沖動頻率都符合分形分布，其函數(shù)關(guān)系為冪律，其形式見公式 (1)。在雙對數(shù)坐標(biāo)下，冪律分布表現(xiàn)為1條斜率為負(fù)數(shù)的直線，這一線性關(guān)系是判斷隨機(jī)變量是否滿足冪律的依據(jù)。

p(x)=Cx-D

(1)

式中：p(x)是具有特征線度為x的物體的數(shù)目；C為比例常數(shù)；D為標(biāo)度指數(shù)也稱分形維數(shù)或分維，是描述分形對象的空間分布格局或內(nèi)在的等級層次等最基本的特征量。

采用分形分布理論擬合流域2009—2013年562次降雨事件70次產(chǎn)沙事件中的降雨產(chǎn)沙累積頻度 (Y) 和降雨量 (L) 間的冪律關(guān)系。

2.3 前期有效降雨量計算

引發(fā)流域產(chǎn)沙的累積降雨量閾值為產(chǎn)沙累積降雨閾值 (Lin)，對于觸發(fā)產(chǎn)沙的降雨，當(dāng)流域累積降雨量達(dá)到Lin值，就有某一產(chǎn)沙累積頻度的可能在該流域內(nèi)引發(fā)產(chǎn)沙，將該產(chǎn)沙累積頻度定義為Ycf。即當(dāng)產(chǎn)沙累積頻度為Ycf時 ，產(chǎn)沙前n天的累積降雨閾值減去前n-1天的累積降雨閾值即為第n天的降雨經(jīng)衰減后對n-1天產(chǎn)沙發(fā)生影響的有效降雨量，產(chǎn)沙發(fā)生前多天降雨衰減后的累計有效降雨量即為流域產(chǎn)沙發(fā)生時的前期有效降雨量，據(jù)此可得流域產(chǎn)沙發(fā)生時的前期有效降雨量，方法如下：

Pean=Linn-Linn-1

(2)

(3)

式中：PEA為前期有效降雨量；Pean為第n天的降雨經(jīng)衰減后對n-1天產(chǎn)沙發(fā)生影響的有效降雨量，Linn為產(chǎn)沙前n天的累積降雨閾值；Linn-1前n-1天的累積降雨閾值，m為前期降雨天數(shù)。

3 結(jié)果與分析

3.1 累積產(chǎn)沙頻度與累積降雨量分形分析

2009—2013年期間流域內(nèi)最長降雨持續(xù)歷時為11 d，考慮導(dǎo)致產(chǎn)沙的降雨歷時時段內(nèi)至其前10天的降雨情況。以侵蝕產(chǎn)沙當(dāng)天作為前第1天，將每次產(chǎn)沙事件劃分前1、2、4、6、8、11天，6個時間段的產(chǎn)沙累積頻度和降雨量繪制成底數(shù)為10的雙對數(shù)坐標(biāo)圖 (圖1)。

圖1不同時段內(nèi)累積產(chǎn)沙頻度-累積降雨分布圖
Fig.1 Distribution cumulative frequency of sediment yield-cumulative rainfall in different periods

從圖1可以看出，6個時段的產(chǎn)沙累積頻度隨降雨量的變化均出現(xiàn)2個不同的無標(biāo)度區(qū)間，分別按照公式 (1) 擬合得到的2條不同斜率 (標(biāo)度指數(shù)) 的直線。在6個時段的產(chǎn)沙累積頻度與降雨量分形擬合中決定系數(shù) (R2) 最大值為0.999 4，最小值為0.782 2，均值為0.924 7，整體擬合效果較好，即表明產(chǎn)沙累積頻度與降雨量之間遵循分形分布。

對于關(guān)嶺享樂小流域的產(chǎn)沙發(fā)生當(dāng)天至其前2、4、6、8、11天，6個時間段，不同標(biāo)度區(qū)間交點處的Y值為58%～78%，平均為68.88%，整體在70%左右，當(dāng)降雨量超過交點處的L值后，產(chǎn)沙的累積頻度僅為22%～42% (平均31.12%)。該點的降雨量為引發(fā)產(chǎn)沙的累積降雨閾值的上邊界，即在該流域引發(fā)70%左右產(chǎn)沙的累積降雨閾值 (Lin) 分別為31.6 mm (1 d)、37.8 mm (2 d)、41.4 mm (4 d)、51.7 mm (6 d)、69.9 mm (8 d)、77.8 mm (11 d)。

3.2 累積降雨閾值與時段分形分析

將流域各時段70%產(chǎn)沙可能時的累積降雨閾值 (Lin) 作為累積降雨時段天數(shù) (n) 的函數(shù)標(biāo)繪為橫坐標(biāo)以2為底數(shù)，縱坐標(biāo)以10為底數(shù)的雙對數(shù)坐標(biāo)圖 (圖2)。

圖2累積降雨閾值與累積降雨時段的分形關(guān)系
Fig.2 A power-exponential relation between cumulative rainfall threshold and cumulative rainfall time

可以發(fā)現(xiàn)兩者之間的函數(shù)遵循冪律分布：

Lin=29.147n0.3741≤n≤11

(4)

3.3 前期有效降雨分析

將侵蝕產(chǎn)沙當(dāng)天作為前1天，對于前n天由公式 (4) 有

Linn=29.147n0.374(1≤n≤m)

(5)

對于前n-1天有

Linn-1=29.147(n-1)0.374

(6)

聯(lián)合式 (5) 有

(7)

聯(lián)合式 (2)、式 (3)、式 (6)、式 (7) 得到對于前n-1天有效的降雨量Pean：

(8)

(9)

對于Lin有

(10)

式中：Pi為相對產(chǎn)沙當(dāng)天前第i天的降雨量，當(dāng)i取1時，P1為產(chǎn)沙當(dāng)天的降雨量。

前期累積降雨對流域產(chǎn)沙發(fā)生的影響即為前期有效降雨量。對于流域侵蝕產(chǎn)沙發(fā)生時 (前第1天) 的有效前期降雨量PEA，由式 (9) 有

(11)

使用有效前期降雨量模式 (11)，計算出關(guān)嶺享樂小流域2009—2013年562個降雨事件的累積前10 d前期有效降雨量。其中日降雨量最大值為102.9 mm，最小值為0.1 mm，當(dāng)日前期有效降雨量為分別13 mm、44.89 mm；前期有效降雨最大為117 mm，并有多天為0 mm。不同降雨日降雨情況異質(zhì)性強(qiáng)，降雨量、前期有效降雨量變異系數(shù)分別為1.69、1.38，降雨量的變異性大于前期有效降雨量，流域日降雨量與前期有效降雨量變化不同步，發(fā)生產(chǎn)沙的日降雨量最小值為13 mm，該天前期有效前期降雨量為29.33 mm。研究期內(nèi)共有106 d降雨量大于13 mm，其中38 d未發(fā)生產(chǎn)沙，未產(chǎn)沙降雨日最大降雨量為47.4 mm，其前期有效降雨為16.65 mm?？芍?，日降雨量或前期有效降雨量單一因素不能預(yù)測流域產(chǎn)沙事件的發(fā)生。

3.4 侵蝕降雨產(chǎn)沙預(yù)測

用日降雨量與前期有效降雨量共同確定流域產(chǎn)沙事件的發(fā)生，將小流域2009—2012年中共452個降雨日的流域日降雨量作為前期有效降雨量的函數(shù)標(biāo)繪在雙對數(shù)坐標(biāo)圖上 (圖3)。連接產(chǎn)沙時的最小日降雨量與最小前期有效降雨點，獲得引發(fā)該流域70%以上的產(chǎn)沙可能的降雨閾值下邊界線 (圖4中的直線)，線性方程為：

P1=20.879-0.282 9PEA

(12)

圖3引發(fā)產(chǎn)沙的最小降雨閾值
Fig.3 Minimum rainfall threshold triggers sediment yield

將閾值下邊界線作為流域產(chǎn)沙降雨條件的判別線，當(dāng)前期有效降雨量與當(dāng)日降雨量點 (PEA,P1) 落在閾值線以上時，視為流域發(fā)生侵蝕性產(chǎn)沙，當(dāng)點 (PEA,P1) 落在閾值邊界線以下時流域不發(fā)生侵蝕產(chǎn)沙。代表降雨日的452個點中有128個落于閾值邊界線上，324個落于閾值邊界線下，與流域?qū)嶋H59 d產(chǎn)沙、303 d未發(fā)生產(chǎn)沙相比，準(zhǔn)確率為86.95%。則有當(dāng)前期有效降雨量與當(dāng)日降雨量點 (PEA,P1) 落在閾值邊界線以上時，該流域發(fā)生侵蝕性產(chǎn)沙可能性較高，當(dāng)點 (PEA,P1) 落在閾值線以下時發(fā)生侵蝕產(chǎn)沙的可能性較低。

將關(guān)嶺享樂小流域2013年發(fā)生的降雨前10天前期有效降雨量與日降雨量標(biāo)繪在雙對數(shù)坐標(biāo)圖上，并添加日降雨閾值下邊界線 (圖4)。最小降雨閾值以上和以下的 (PEA,P1) 點分別為9個與101個，即2013年該流域應(yīng)有9次侵蝕產(chǎn)沙降雨與101次非侵蝕性降雨。對比發(fā)現(xiàn)該流域2013年實際共有110 d發(fā)生降雨，其中侵蝕產(chǎn)沙型降雨11次，非侵蝕型降雨99次，與預(yù)測值相差4次。在70%侵蝕產(chǎn)沙可能情況下，預(yù)測準(zhǔn)確106次，預(yù)測準(zhǔn)確率達(dá)到96.36%。對于享樂小流域，將日降雨閾值下邊界線與產(chǎn)沙發(fā)生有關(guān)的區(qū)域平均降雨條件聯(lián)合能很好的預(yù)測流域的侵蝕產(chǎn)沙的發(fā)生。

圖42013年降雨-產(chǎn)沙情況
Fig.4 The relationship of rainfall and sediment in 2013

4 結(jié)論與討論

1) 關(guān)嶺享樂小流域地表侵蝕產(chǎn)沙事件的累積頻度和降雨量、產(chǎn)沙事件的累積降雨閾值和降雨時段之間均遵循分形分布的冪指數(shù)關(guān)系。利用分形分布理論結(jié)合流域降雨、產(chǎn)沙情況可建立適用于喀斯特地區(qū)的前期有效降雨模型。

2) 基于累積產(chǎn)沙頻度-降雨量分形關(guān)系建立的前期有效降雨模型考慮了降雨衰減系數(shù)并非經(jīng)驗系數(shù)，而是由流域內(nèi)引發(fā)產(chǎn)沙的累積降雨閾值隨觀測時段天數(shù)變化的標(biāo)度指數(shù)所決定，能去除衰減系數(shù)的取值隨意性，相比傳統(tǒng)模型可較為符合實際的計算出產(chǎn)沙發(fā)生時的流域前期有效降雨量。

3) 將日降雨量與有效前期降雨量結(jié)合可以確定喀斯特流域發(fā)生侵蝕性降雨產(chǎn)沙的降雨閾值。關(guān)嶺享樂小流域基于前期有效降雨量PEA、當(dāng)日降雨量P1的流域地表侵蝕產(chǎn)沙發(fā)生閾值線性方程為P1=20.879-0.282 9PEA，當(dāng)前期有效降雨量與當(dāng)日降雨量點 (PEA,P1) 落在閾值線以下不發(fā)生流域地表侵蝕產(chǎn)沙，落在閾值范圍以上發(fā)生侵蝕產(chǎn)沙。該閾值線性方程對關(guān)嶺享樂小流域2013年產(chǎn)沙事件的預(yù)測準(zhǔn)確率達(dá)到96.36%。

4) 對喀斯特流域，通過降雨閾值根據(jù)日降雨預(yù)報即可得到有效、較為準(zhǔn)確的流域產(chǎn)沙預(yù)測結(jié)果。因此，可將前期有效降雨量作為是否發(fā)生產(chǎn)沙的一個重要參數(shù)，對喀斯特流域侵蝕產(chǎn)沙發(fā)生的預(yù)測具有重要意義。

目前前期有效降雨計算方法多是是基于蒸發(fā)量、降雨量等基本指標(biāo)構(gòu)建的簡單模型[4-7]，但是各自所取的衰減系數(shù)因各地區(qū)太熱輻射、氣溫等小生境氣候因子不同而具有較大差異，多采取區(qū)域經(jīng)驗值，具有隨意性。采用累積產(chǎn)沙頻度-降雨量分形關(guān)系計算該流域前期有效降雨量，能去除衰減系數(shù)的取值隨意性。但因各流域的產(chǎn)沙情況與日降雨情況不同，其標(biāo)度指數(shù)也不同，對于不同流域需各自計算其標(biāo)度指數(shù)，構(gòu)建前期有效降雨模型，相比傳統(tǒng)模型計算量大且復(fù)雜，對該方法的應(yīng)用有一定的難度。通過建立多個喀斯特流域的累積產(chǎn)沙頻度-降雨量前期有效降雨模型計算更為精確、適用于喀斯特地區(qū)的經(jīng)驗衰減系數(shù)可進(jìn)一步的研究。

基于分形分布計算自然災(zāi)害預(yù)測被較多的應(yīng)用于常態(tài)地貌的滑坡、泥石流等地質(zhì)災(zāi)害預(yù)測中，對于喀斯特地區(qū)鮮見?？λ固氐貐^(qū)由于復(fù)雜的地形地貌和地下、地表二元空間結(jié)構(gòu)影響，各種環(huán)境因子之間存在復(fù)雜的相互關(guān)系，釆用常規(guī)分析方法很難準(zhǔn)確界定各因子的影響程度，難以預(yù)測土壤侵蝕的發(fā)生，更難以有效預(yù)測滑坡、崩塌、泥石流、巖溶塌陷等自然災(zāi)害的發(fā)生。分形分布是普遍存在于自然界的，所有規(guī)模的沖動頻率都符合分形分布，本研究通過利用分形分布計算出流域產(chǎn)沙閾值，對流域侵蝕產(chǎn)沙做出很好的參考，基于分形分布計算自然災(zāi)害發(fā)生前環(huán)境因子狀態(tài)，也可以建立災(zāi)害發(fā)生時的影響因子自然臨界閾值，對今后在喀斯特地區(qū)自然災(zāi)害的預(yù)測提供參考。

[參 考 文 獻(xiàn)]

[1] 凡非得, 王克林, 熊鷹, 等. 西南喀斯特區(qū)域水土流失敏感性評價及其空間分異特征[J]. 生態(tài)學(xué)報, 2011, 31(21): 6353-6362.

[2] 曹建華, 袁道先, 章程, 等. 受地質(zhì)條件制約的中國西南巖溶生態(tài)系統(tǒng)[J]. 地球與環(huán)境, 2004, 32(1): 1-8.

[3] 李瑞, 李勇, 劉云芳. 貴州喀斯特地區(qū)降雨與坡面土壤侵蝕關(guān)系研究[J]. 水土保持研究, 2012, 19(3): 7-11.

[4] 劉正堂, 戴全厚, 倪九派, 等. 喀斯特地區(qū)裸坡面土壤侵蝕的人工模擬降雨試驗研究[J]. 水土保持學(xué)報, 2013, 27(5): 12-16.

[5] Glade T, Crazier M, Smith P. Applying probability determination to refine landslide-triggering rainfall thresholds using an empirical “Antecedent Daily Rainfall Model” [J]. Pure and Applied Geophysics, 2000, 157(6): 1059-1079.

[6] Zêzere J L, Trigo R M, Trigo I F. Shallow and deep landslides induced by rainfall in the Lisbon region (Portugal): assessment of relationships with the North Atlantic Oscillation[J]. Natural Hazards and Earth System Science, 2005, 5(3): 331-344.

[7] 姜志文, 劉樹春. 柴河水庫流域前期影響雨量的計算方法[J]. 現(xiàn)代農(nóng)業(yè)科技, 2010(9): 258, 260.

[8] 崔鵬, 楊坤, 陳杰. 前期降雨對泥石流形成的貢獻(xiàn): 以蔣家溝泥石流形成為例[J]. 中國水土保持科學(xué), 2003, 1(1): 11-15.

[9] 李昌志, 劉興年, 曹叔尤, 等. 前期降雨與不同沙源條件小流域產(chǎn)沙關(guān)系的對比研究[J]. 水土保持學(xué)報, 2001, 15(2): 36-39.

[10] 陳洪松, 王克林. 西南喀斯特山區(qū)土壤水分研究[J]. 農(nóng)業(yè)現(xiàn)代化研究, 2008, 29(6): 734-738.

[11] 沈紫燕, 王輝, 平李娜, 等. 前期土壤含水量對粘性紅壤產(chǎn)流產(chǎn)沙及溶質(zhì)運移的影響[J]. 水土保持學(xué)報, 2014, 28(1): 58-62.

[12] 王海斌. 不同水土保持措施下徑流小區(qū)降雨與產(chǎn)流產(chǎn)沙關(guān)系研究[J]. 水土保持研究, 2011, 18(5): 63-66.

[13] 紀(jì)啟芳, 張興奇, 張科利, 等. 貴州省喀斯特地區(qū)坡面產(chǎn)流產(chǎn)沙特征[J]. 水土保持研究, 2012, 19(4): 1-5.

[14] 李長江, 麻土華, 孫樂玲, 等. 降雨型滑坡預(yù)報中計算前期有效降雨量的一種新方法[J]. 山地學(xué)報, 2011, 29(1): 81-86.

[15] 李長江, 麻土華, 朱興盛. 降雨型滑坡預(yù)報的理論、方法及應(yīng)用[M]. 北京: 地質(zhì)出版社, 2008.