本發(fā)明涉及一種水土保持工程減水減沙效益定量評價方法，是一種水文計算方法，是一種使用計算機(jī)技術(shù)的水文模擬和計算方法。

背景技術(shù)：

現(xiàn)有的水土保持水文水資源效應(yīng)研究所采取的方法均為“水文法”和“水保法”。由于這 兩種方法都不能系統(tǒng)地揭示水土保持對水循環(huán)影響的物理機(jī)制,并且采用的水資源評價口徑。只是單一的狹義水資源,使得研究成果不能真實反映水土保持的水文水資源效應(yīng)。以狹義水資源和廣義水資源為評價口徑的基于物理機(jī)制的分布式 水文模型則能克服水文法和水保法的嚴(yán)重弊端。

“水文法”也叫統(tǒng)計分析法，利用流域水文泥沙觀測資料分析水土保持措施蓄水?dāng)r沙作用的一種方法。基本原理是基于對降雨產(chǎn)流產(chǎn)沙規(guī)律的分析，建立水土保持治理前流域產(chǎn)流產(chǎn)沙關(guān)系式，然后將治理后的降雨條件帶入關(guān)系式，計算得到該降雨條件時流域未治理條件下的產(chǎn)水產(chǎn)沙量，再與治理后的實測產(chǎn)水產(chǎn)沙量比較，從而求得水土保持措施對流域水沙的影響量。

“水保法”也叫成因分析法，是根據(jù)不同水土保持措施減水減沙作用實測成果，考慮產(chǎn)沙在河道中的輸移及沖淤變化，以及新增水土流失數(shù)量等，計算水土保持減水減沙效益的一種方法。一般直接采用具體水土保持措施的面積和該項措施的蓄水保土定額之積作為該項水土保持措施的減水減沙效益。

上述兩種方法的問題都在于不能精確的計算或模擬坡面的水沙過程，因此，大大的影響了計算或模擬的精確度。

技術(shù)實現(xiàn)要素：

為了克服現(xiàn)有技術(shù)的問題，本發(fā)明提出了一種水土保持工程減水減沙效益定量評價方法。所述的方法基于坡面或流域水沙物理過程，或者將不同分辨率的柵格作為計算單元，或者將柵格整合為等高帶作為計算單元，建立計算單元條件下的降雨-產(chǎn)流-匯流和土壤侵蝕-輸移機(jī)制的水沙過程，用以水土保持工程減水減沙效益定量評價的模擬方法。

本發(fā)明的目的是這樣實現(xiàn)的：一種水土保持工程減水減沙效益定量評價方法，所述方法的步驟如下：

地形數(shù)字化處理的步驟：用于DEM數(shù)據(jù)為基礎(chǔ)，以ArcGIS水文模塊為主要工具，采用D8法，建立以柵格或等高帶為基本的計算單元的坡面水文過程計算所需的地形數(shù)據(jù)、匯流關(guān)系數(shù)據(jù)文件，以這些文件為基礎(chǔ)，建立坡面侵蝕與泥沙輸移過程計算所需的地形數(shù)據(jù)文件；

水土保持措施數(shù)據(jù)收集的步驟：用于收集模擬區(qū)域內(nèi)不同時間的坡面植被、土壤、土地利用情況，以及水土保持措施情況，形成不同時間段的水土保持?jǐn)?shù)據(jù)集；

其他數(shù)據(jù)收集與處理的步驟：用于收集模擬區(qū)域內(nèi)及其附近的水文、氣象，其中水文氣象數(shù)據(jù)包括降雨、平均風(fēng)速、平均氣溫、日照時數(shù)、空氣相對濕度，對收集到的數(shù)據(jù)按照模擬計算所需的時間步長進(jìn)行處理和空間匹配處理；

水沙計算的步驟：用于輸入模型率定參數(shù)或應(yīng)用參數(shù)，進(jìn)行坡面水文過程計算、坡面侵蝕與泥沙輸移過程計算：

所述的坡面水文過程計算為：從坡面頂端的計算單元開始，進(jìn)行計算單元中坡面水文過程計算，包括：蒸發(fā)蒸騰、入滲、地表徑流、壤中徑流、坡面匯流、積雪融雪計算，如果存在上游計算單元則接收上游計算單元的參數(shù)，加入到計算中；

所述的坡面侵蝕與泥沙輸移過程計算為：從坡面頂端的計算單元開始，進(jìn)行計算單元中坡面侵蝕與泥沙輸送和移動過程模擬和計算，所述的模擬和計算包括：雨滴濺蝕模擬、薄層水流侵蝕模擬、股流侵蝕過程模擬、重力侵蝕過程模擬，如果存在上游計算單元則接收上游計算單元的參數(shù)，加入到計算中；

判斷的步驟：根據(jù)匯流關(guān)系文件給定的計算關(guān)系， 進(jìn)行過程計算和決定是否輸出當(dāng)前計算結(jié)果，如果“是”進(jìn)入“計算偏差和率定參數(shù)的步驟”，如果“否”則進(jìn)入“傳遞變量參數(shù)的步驟”；

模型參數(shù)的率定和模型驗證的步驟：如果輸入的參數(shù)為模型率定參數(shù)，則通過計算偏差對水沙計算的準(zhǔn)確度進(jìn)行評價，以率定模型是否符合真實情況，并不斷的調(diào)整參數(shù)并回到“水沙計算的步驟”進(jìn)行計算，同時不斷的進(jìn)行偏差計算和比較，使模型計算與實際情況的偏差達(dá)到最小，并進(jìn)行模型的驗證，如果輸入的參數(shù)為應(yīng)用參數(shù)則跳過本步驟；

傳遞變量參數(shù)的步驟：用于收集“水沙計算的步驟”所獲得的變量參數(shù)，并攜帶這些變量參數(shù)回到“水沙計算的步驟”，繼續(xù)進(jìn)行水沙計算；

設(shè)置水土保持情景的步驟：用于使用水土保持措施數(shù)據(jù)，設(shè)置多個年代的水土保持措施情景；

生成不同水土保持情景下的水沙過程數(shù)據(jù)的步驟：用于將設(shè)計的水土保持情景數(shù)據(jù)運(yùn)行計算，得到多個水土保持措施的水沙過程數(shù)據(jù)；

分析計算不同水土保持情景的保水保土效益的步驟：用于對比不同水土保持情景下的水沙過程變化量進(jìn)行分析；

比較并結(jié)束的步驟：用于對各個情景水土保持情景的水土保持?jǐn)?shù)據(jù)集所對應(yīng)的水沙過程時間段數(shù)據(jù)進(jìn)行比較，并輸出比較結(jié)果，結(jié)束計算過程。

進(jìn)一步的，所述的蒸發(fā)蒸騰計算：

式中，F_W、F_U、F_SV、F_IR、F_NI分別為計算單元內(nèi)水域、不透水域、裸地-植被域、灌溉農(nóng)田及非灌溉農(nóng)田的面積率；E_W、E_SV、E_U、E_IR、E_NI分別為計算單元內(nèi)水域、不透水域、裸地－植被域、灌溉農(nóng)田及非灌溉農(nóng)田的蒸發(fā)量或蒸發(fā)蒸騰量。

進(jìn)一步的，所述的入滲計算采用Green-Ampt鉛直一維入滲模型模擬降雨入滲及超滲坡面徑流，和通用Green-Ampt模型進(jìn)行計算。

進(jìn)一步的，所述的地表徑流和壤中徑流的計算，

地表徑流為：地表徑流等于降雨減去降雨時的蒸發(fā)損失，

不透水域的地表徑流按上述公式所述的地表徑流和壤中徑流的計算，

地表徑流為：地表徑流等于降雨減去降雨時的蒸發(fā)損失，

不透水域的地表徑流：；

及

計算，

式中，P為降雨，H_u為洼地儲蓄，E_u為蒸發(fā)，R_u為表面徑流，H_umax為最大洼地儲蓄深，E_umax為潛在蒸發(fā)，c為城市建筑物在不透水域的面積率，下標(biāo)1表示城市建筑物，下標(biāo)2表示城市地表面；

裸地－植被域的地表徑流則根據(jù)降雨強(qiáng)度是否超過土壤的入滲能力分以下兩種情況計算：霍頓坡面徑流和飽和坡面徑流，

壤中徑流R2由下式計算：

式中，k(θ)為體積含水率q對應(yīng)的沿山坡方向的土壤導(dǎo)水系數(shù)，slope為地表面坡度，L為計算單元內(nèi)的河道長度，d為不飽和土壤層的厚度。

進(jìn)一步的，所述的坡面匯流計算為：

運(yùn)動波方程：

式中，A為流水?dāng)嗝婷娣e，Q為斷面流量，q_L為網(wǎng)格單元或河道的單寬流入量，n為曼寧糙率系數(shù)，R為水力半徑，S₀為網(wǎng)格單元地表面坡降或河道的縱向坡降，S_f為摩擦坡降，

動力波方程：

式中，V為斷面流速，V_x為單寬流入量的流速在x方向的分量。

進(jìn)一步的，所述的積雪融雪過程計算：

式中，SM為融雪量，M_f為融化系數(shù)，T_a為氣溫指標(biāo)，T₀為融化臨界溫度，S為積雪水當(dāng)量，SW為降雪水當(dāng)量，E_snow為積雪升華量。

進(jìn)一步的，所述的雨滴濺蝕模擬和計算：

式中，D₁為雨滴擊濺侵蝕量，E_drop為雨滴動能，I為雨強(qiáng)，J₁為地表坡度，k₁，α₁，β₁為經(jīng)驗參數(shù)。

其中雨動能E的計算：

式中，E_unit為單位降雨動能，k₁＇,α₁＇為經(jīng)驗參數(shù)；

當(dāng)水深大于雨滴直徑3倍以上時，取水深大于0.6cm時，雨滴濺蝕作用消失；

由雨滴濺蝕增加的土壤侵蝕輸沙能力計算公式為：

式中， q_s1為單寬流量輸沙能力，k₂，α₂，β₂為經(jīng)驗參數(shù)。

進(jìn)一步的，所述的薄層水流侵蝕模擬：

式中， Dc為水流剝離土壤速率，k₃為土壤的可蝕性參數(shù)，τ_f為水流對土壤顆粒的剪切應(yīng)力，τ_c為土壤的臨界抗剪切應(yīng)力，Dr為細(xì)溝水流剝蝕率，q為單寬流量，c為水流泥沙含量， Tc為水流的挾沙能力，k₄，α₄為經(jīng)驗常數(shù)。

進(jìn)一步的，所述的股流侵蝕過程模擬，股流的挾沙能力T_SE計算公式為：

，

式中，k₅為淺溝水流攜沙能力系數(shù)，m為側(cè)向匯流影響常數(shù)，ω_u為單位水流功率；

式中，S_E為股流侵蝕量， Q_E為流量， D_r為上游來沙量。

進(jìn)一步的，所述的重力侵蝕過程模擬，重力侵蝕量V_g計算公式為：

式中，k₆為發(fā)生重力侵蝕的溝長系數(shù)，L_gully為溝道總長度。

本發(fā)明產(chǎn)生的有益效果是：本發(fā)明通過收集坡面的水土保持措施的數(shù)據(jù)以及其他水文、氣象數(shù)據(jù)，利用不同分辨率的柵格DEM數(shù)據(jù)或等高帶數(shù)據(jù)為平臺，以山坡真實水沙物理過程為圖景，實現(xiàn)了坡面降雨-徑流和侵蝕-泥沙過程的耦合模擬計算能力，在精確水沙過程計算的前提下，對不同時期的坡面水土保持狀況進(jìn)行比較分析，十分精確的獲得水土保持措施的效果。本發(fā)明由于精確的模擬和計算水沙過程與已有技術(shù)相比，對坡面的水土保持狀況的分析更加精確，更加接近實際狀況，對如何設(shè)計水土保持措施提供了更加切合實際的科學(xué)依據(jù)。

附圖說明

下面結(jié)合附圖和實施例對本發(fā)明作進(jìn)一步說明。

圖1是本發(fā)明的實施例一所述方法的流程示意圖；

圖2是本發(fā)明的實施例一所述應(yīng)用實例中涇河流域降雨-徑流量年代變化；

圖3是本發(fā)明的實施例一所述應(yīng)用實例中涇河流域降雨-輸沙量年代變化；

圖4是本發(fā)明的實施例一所述應(yīng)用實例中不同情景下流域單位徑流輸沙量；

圖5是本發(fā)明的實施例十所述重力侵蝕原理示意圖。

具體實施方式

實施例一：

本實施例是一種水土保持工程減水減沙效益定量評價方法。本實施例以基于過程的水文模擬為基礎(chǔ)，通過系統(tǒng)識別地形特征與水流動力學(xué)特征的耦合關(guān)系，建立以地形“面（片）蝕-細(xì)溝侵蝕-淺溝侵蝕-切溝侵蝕”土壤侵蝕鏈為平臺，以“薄層水流-股流”為典型水動力條件的土壤侵蝕與泥沙輸移過程模擬方法。本實施例首先通過對DEM數(shù)據(jù)的處理形成坡面計算單元拓?fù)潢P(guān)系和匯流溝道網(wǎng)，然后此為平臺對水文氣象、水土保持措施（水土保持林草、梯田、水平溝、魚鱗坑）、土地利用、地形等信息進(jìn)行處理，將處理后的數(shù)據(jù)輸入模型，接著進(jìn)行產(chǎn)匯流計算和侵蝕產(chǎn)沙輸沙計算，最后利用已有的水文泥沙數(shù)據(jù)對模型進(jìn)行率定和校核。從而實現(xiàn)坡面水沙過程模擬。所述實施例綜合構(gòu)建侵蝕地貌形態(tài)和水動力學(xué)過程特點的流域分布式水沙耦合模型。針對坡面股流侵蝕和重力侵蝕研究相對不足的現(xiàn)狀，通過室內(nèi)試驗資料分析，以及黃土抗剪強(qiáng)度變化規(guī)律野外試驗對黃土區(qū)股流侵蝕過程輸沙能力以及重力侵蝕中土壤抗剪強(qiáng)度變化規(guī)律兩個基本問題進(jìn)行了研究。從而形成了物理機(jī)制相對完善的分布式流域水沙耦合模型。

所述方法的步驟如下（流程示意圖見圖1）：

一、地形數(shù)字化處理的步驟：用于DEM數(shù)據(jù)為基礎(chǔ)，以ArcGIS水文模塊為主要工具，采用D8法，建立以柵格為基本計算單元的坡面水文過程計算所需的地形數(shù)據(jù)、匯流關(guān)系數(shù)據(jù)文件，以這些文件為基礎(chǔ)，建立坡面侵蝕與泥沙輸移過程計算所需的地形數(shù)據(jù)文件。

（1）以高分辨率柵格DEM數(shù)據(jù)為基礎(chǔ)，以ArcGIS水文模塊為主要工具，采用D8法，建立以柵格為基本計算單元的坡面水文過程計算所需的地形數(shù)據(jù)、匯流關(guān)系數(shù)據(jù)文件。

（2）以上述文件為基礎(chǔ)，按照下述方法建立坡面侵蝕與泥沙輸移過程計算所需的地形數(shù)據(jù)文件。

不同侵蝕地貌的坡面分布密度運(yùn)用不同侵蝕輸沙規(guī)律的基礎(chǔ)。細(xì)溝是黃土坡面分布最廣的溝蝕類型之一。由于細(xì)溝密度和深度在坡面的分布隨坡面長度的增大呈現(xiàn)多峰變化的規(guī)律，分布密度最大可以達(dá)到6%，深度變化在0-14cm之間?？梢岳貌煌恋乩脳l件下細(xì)溝密度概念進(jìn)行等高帶內(nèi)淺溝侵蝕計算。具體模擬過程中假定面（片）蝕水流在整個計算單元內(nèi)發(fā)生，計算單元內(nèi)依據(jù)土地利用的不同，分別給出相對于裸地的衰減系數(shù)；細(xì)溝侵蝕過程模擬中細(xì)溝尺寸依據(jù)不同土地利用類型按面積比進(jìn)行概化。淺溝和切溝根據(jù)具體的地形條件決定是否發(fā)生侵蝕及溝道數(shù)量。概化的典型侵蝕地貌單元的斷面參數(shù)及其發(fā)生的臨界地形和水動力學(xué)條件。

二、水土保持措施數(shù)據(jù)收集的步驟：用于收集模擬區(qū)域內(nèi)不同時間的坡面植被、土壤、土地利用情況，以及水土保持措施情況，形成不同時間段的水土保持?jǐn)?shù)據(jù)集。

三、其他數(shù)據(jù)收集與處理的步驟：其他數(shù)據(jù)收集與處理的步驟：用于收集模擬區(qū)域內(nèi)及其附近的水文、氣象，其中水文氣象數(shù)據(jù)包括降雨、平均風(fēng)速、平均氣溫、日照時數(shù)、空氣相對濕度，對收集到的數(shù)據(jù)按照模擬計算所需的時間步長進(jìn)行處理和空間匹配處理，同時進(jìn)入“坡面水文過程計算的步驟”和“坡面侵蝕與泥沙輸移過程計算的步驟”。

數(shù)據(jù)收集：收集模擬區(qū)域內(nèi)及其附近的水文、氣象、植被、土壤、土地利用、水土保持措施數(shù)據(jù)，其中水文氣象數(shù)據(jù)包括降雨、平均風(fēng)速、平均氣溫、日照時數(shù)、空氣相對濕度；植被數(shù)據(jù)包括植被覆蓋度、葉面積指數(shù)；土壤數(shù)據(jù)包括土壤類型、土壤厚度；土地利用數(shù)據(jù)包括模擬時段內(nèi)多期土地利用數(shù)據(jù)，一般5年需要一期；水土保持?jǐn)?shù)據(jù)包括魚鱗坑、水平溝、梯田、水土保持林草、淤地壩等水土保持措施實施的時間、地點、范圍。

數(shù)據(jù)處理：

（1）時間過程處理：對上述收集到的數(shù)據(jù)按照模擬計算所需的時間步長進(jìn)行處理。對于不同數(shù)據(jù)和具體的時間過程有不同的處理方法，相關(guān)文獻(xiàn)非常多。

（2）空間匹配處理：首先將上述數(shù)據(jù)進(jìn)行空間差值和展布，形成覆蓋模擬區(qū)域的空間分布數(shù)據(jù)，然后與數(shù)字化坡面的計算單元即柵格進(jìn)行匹配。相關(guān)文獻(xiàn)非常多。對于坡面侵蝕過程，根據(jù)匯流條件將計算單元劃分為由平面、細(xì)溝、淺溝和切溝構(gòu)成的亞計算單元地貌形態(tài)。根據(jù)已有的野外調(diào)查數(shù)據(jù)，首先對細(xì)溝、淺溝和切溝的形態(tài)進(jìn)行研究，概括出不同侵蝕地貌的“典型侵蝕形態(tài)”概念。

以上步驟是針對整個模擬區(qū)域進(jìn)行的地形的數(shù)字化和數(shù)據(jù)的處理。以下步驟則是針對一個計算單元所進(jìn)行的計算和模擬。當(dāng)一個計算單元計算完成后則對該計算單元下游的計算單元進(jìn)行計算和模擬，從坡頂一直計算到坡底。一個坡面有多個這樣從坡頂?shù)狡碌椎挠嬎銌卧蛄?，可以同時進(jìn)行各個計算單元序列的計算，也可以在計算完一個計算單元序列后在計算相鄰的計算單元序列，以便應(yīng)用兩個相鄰計算單元需要的變量參數(shù)傳遞。

四、水沙計算的步驟：用于輸入模型率定參數(shù)或應(yīng)用參數(shù)，進(jìn)行坡面水文過程計算、坡面侵蝕與泥沙輸移過程計算：

所述的模型率定參數(shù)是指在模型建立過程中，首先要對模型是否符合實際請進(jìn)行率定，因此先要輸入一些結(jié)果已知的已有參數(shù)對模型進(jìn)行調(diào)整和校驗，并不斷的進(jìn)行偏差計算，不斷的獲取模型計算結(jié)果與實際結(jié)果的偏離度，使模型的計算結(jié)果與實際結(jié)果的偏離度最小，以此建立十分精確的模型，這些率定模型的已知結(jié)果的參數(shù)稱之為“模型率定參數(shù)”。

所述的應(yīng)用參數(shù)是在模型經(jīng)過率定后，已經(jīng)到達(dá)了與實際十分接近的程度，則輸入一些預(yù)測或其他應(yīng)用的參數(shù)，以獲取一個未知的結(jié)果，這些參數(shù)稱之為“應(yīng)用參數(shù)”。

所述的坡面水文過程計算為：從坡面頂端的計算單元開始，進(jìn)行計算單元中坡面水文過程計算，包括：蒸發(fā)蒸騰、入滲、地表徑流、壤中徑流、坡面匯流、積雪融雪計算，如果存在上游計算單元則接收上游計算單元的參數(shù)，加入到計算中。

坡面水沙過程的模擬是基于物理機(jī)制分布式水沙耦合模型的基礎(chǔ)。從侵蝕動力學(xué)角度，黃土區(qū)坡面典型水沙過程主要包括雨滴濺蝕過程、薄層水流侵蝕過程、股流侵蝕過程和重力侵蝕過程。其中面（片）蝕和細(xì)溝侵蝕適用于薄層水流侵蝕模擬，淺溝侵蝕和切溝侵蝕適用于股流侵蝕。目前對雨滴濺蝕和薄層水流侵蝕過程的模擬研究已經(jīng)較為成熟，而對股流侵蝕和重力侵蝕過程的機(jī)理描述尚不完善。本實施例采用雨滴濺蝕和薄層水流侵蝕過程模，同時為了建立完善的坡面水沙過程模擬模型，對股流侵蝕過程和重力侵蝕過程在物理圖景概化的基礎(chǔ)上，通過實驗研究建立了基于物理機(jī)制的模擬子模型。

蒸發(fā)蒸騰的計算：可以采用“馬賽克”結(jié)構(gòu)考慮網(wǎng)格單元內(nèi)的土地利用變異問題，每個網(wǎng)格單元的蒸發(fā)蒸騰可能包括植被截留蒸發(fā)、土壤蒸發(fā)、水面蒸發(fā)和植被蒸騰等多項。參照土壤-植被-大氣通量交換方法(SVATS)中的ISBA模型，采用Penman公式或Penman-Monteith公式等進(jìn)行計算。同時，由于蒸發(fā)蒸騰過程和能量交換過程客觀上融為一體，為計算蒸發(fā)蒸騰，地表附近的輻射、潛熱、顯熱和熱傳導(dǎo)的計算不可缺少，而這些熱通量又均是地表溫度的函數(shù)。為減輕計算負(fù)擔(dān)，模型對熱傳導(dǎo)及地表溫度的計算采用強(qiáng)制復(fù)原法。

入滲的計算：可以采用Green-Ampt鉛直一維入滲模型模擬降雨入滲及超滲坡面徑流。同時，考慮到由自然力和人類活動(如農(nóng)業(yè)耕作)等引起的土壤分層問題，模型采用Jia 和Tamai提出的實際降雨條件下的多層Green-Ampt模型，以下稱通用Green-Ampt模型。

地表徑流的計算：水域的地表徑流等于降雨減去降雨時的蒸發(fā)損失，裸地－植被域（透水域）的地表徑流則根據(jù)降雨強(qiáng)度是否超過土壤的入滲能力分以霍頓坡面徑流飽和坡面徑流進(jìn)行計算。

壤中徑流計算：在山地丘陵等地形起伏地區(qū)，同時考慮坡向壤中徑流及土壤滲透系數(shù)的各向變異性。

坡面匯流計算：在河網(wǎng)水系生成和流域劃分的基礎(chǔ)上，根據(jù)網(wǎng)格單元DEM及土地利用等基本信息準(zhǔn)備各子流域等高帶的屬性表（包括面積、長度、寬度、平均高程、坡度和曼寧糙率等），采用一維運(yùn)動波法從上游等高帶至下游等高帶計算坡面匯流，并將下游等高帶的坡面匯流輸入給所在子流域內(nèi)的河道。

積雪融雪計算：可以采用“溫度指標(biāo)法”計算積雪融雪的日變化過程，或其他類似的方法。

所述的坡面侵蝕與泥沙輸移過程計算為：按照匯流關(guān)系數(shù)據(jù)文件從坡面頂端的計算單元開始，進(jìn)行計算單元中坡面侵蝕與泥沙輸送和移動過程模擬和計算，所述的模擬和計算包括：雨滴濺蝕模擬、薄層水流侵蝕模擬、股流侵蝕過程模擬、重力侵蝕過程模擬，如果存在上游計算單元則接收上游計算單元的參數(shù)，加入到計算中。

坡面水沙過程的模擬是基于物理機(jī)制分布式水沙耦合模型的基礎(chǔ)。從侵蝕動力學(xué)角度，黃土區(qū)坡面典型水沙過程主要包括雨滴濺蝕過程、薄層水流侵蝕過程、股流侵蝕過程和重力侵蝕過程。其中面（片）蝕和細(xì)溝侵蝕適用于薄層水流侵蝕模擬，淺溝侵蝕和切溝侵蝕適用于股流侵蝕。目前對雨滴濺蝕和薄層水流侵蝕過程的模擬研究已經(jīng)較為成熟，而對股流侵蝕和重力侵蝕過程的機(jī)理描述尚不完善。

雨滴濺蝕模擬：對于土壤性質(zhì)基本一致的流域，影響降雨濺蝕的主要因素為降雨動能、地表坡度。柵格內(nèi)可依據(jù)土地利用的不同，分別給出相對于裸地的衰減系數(shù)。其中對于裸土地雨滴濺蝕的計算選取適合于黃土地區(qū)，可運(yùn)用綜合考慮降雨及坡度的模型—吳普特模型進(jìn)行濺蝕計算。

薄層水流侵蝕模擬：可采用Foster和Meyer年提出的關(guān)系式描述薄層水流土壤分離速率與輸沙率之間的關(guān)系。該關(guān)系式被用于WEPP模型，且經(jīng)實驗驗證這一假定符合黃土區(qū)薄層水流侵蝕過程規(guī)律。

股流侵蝕過程模擬：股流挾沙能力是指水流能量全部用于泥沙輸移時的水流最大挾沙能力。如果使用經(jīng)過篩分除去石礫的松散土壤，土壤的抗蝕能力下降較自然土壤明顯下降，由于黃土的易蝕性，水流的侵蝕耗能很少減少，在這種情況下水流含沙量可以近似看做股流挾沙能力。

重力侵蝕過程模擬：坡面水沙過程中的重力侵蝕，其影響因素相對單一。對于坡面坡度大于黃土休止角的土體，其是否發(fā)生重力侵蝕取決于土體的下滑力G_x與土體抗剪強(qiáng)度τ_c之間的大小關(guān)系。當(dāng)G_x＞τ_c重力侵蝕發(fā)生。

五、判斷的步驟：用于根據(jù)匯流關(guān)系文件給定的計算關(guān)系， 進(jìn)行過程計算和決定是否輸出當(dāng)前計算結(jié)果，如果“是”進(jìn)入“計算偏差和率定參數(shù)的步驟”，如果“否”則進(jìn)入“傳遞變量參數(shù)的步驟”；

由于計算單元的DEM數(shù)據(jù)是從坡面到坡底順序排列，上下游計算單元之間有變量參數(shù)傳遞的關(guān)系，因此，計算開始時，首先應(yīng)當(dāng)從坡頂開始計算，然后順序向下計算，并不斷的將變量參數(shù)向下一個計算單元傳遞，所以一個計算單元的計算完畢后就要判斷一下，是否達(dá)到了流域出口，如果沒有達(dá)到流域出口，就要繼續(xù)進(jìn)行計算。

六、傳遞變量參數(shù)的步驟：用于收集“水沙計算的步驟”和所獲得的變量參數(shù)，并攜帶這些變量參數(shù)回到“水沙計算的步驟”，繼續(xù)進(jìn)行水沙計算。

由于水流在坡面流動具有從上往下并不斷擴(kuò)大流動范圍的特點，水沙過程的計算中具有下游的計算單元對該計算單元的上游計算單元，以及上游計算單元相鄰的計算單元有數(shù)據(jù)承傳的關(guān)系，因此，當(dāng)一個計算單元完成水沙過程的模擬后，通常需要將水沙過程的計算結(jié)果向下游計算單元傳遞，除非達(dá)到坡底的匯流。

七、模型參數(shù)的率定和模型驗證的步驟：如果輸入的參數(shù)為模型率定參數(shù)，則通過計算偏差對水沙計算的準(zhǔn)確度進(jìn)行評價，以率定模型是否符合真實情況，并不斷的調(diào)整參數(shù)并回到“水沙計算的步驟”進(jìn)行計算，同時不斷的進(jìn)行偏差計算和比較，使模型計算與實際情況的偏差達(dá)到最小，并進(jìn)行模型的驗證，如果輸入的參數(shù)為應(yīng)用參數(shù)則跳過本步驟。

模型參數(shù)的率定和模型驗證是一個反復(fù)進(jìn)行的過程，需要不斷的進(jìn)行偏差計算和調(diào)整、修正參數(shù)。參數(shù)的調(diào)整十分復(fù)雜，需要調(diào)整的參數(shù)往往很多，調(diào)整哪些些參數(shù)和調(diào)整的幅度控制多少，才能盡快準(zhǔn)確的達(dá)到目的，需要智慧和耐力。

模型參數(shù)率定過程是按照模型框架需求，對模型參數(shù)進(jìn)行調(diào)整使模型準(zhǔn)確反映研究對象真實世界客觀規(guī)律的過程。模型驗證則是對模型參數(shù)率定結(jié)果的可靠性和準(zhǔn)確性檢驗。驗證標(biāo)準(zhǔn)一般有Nash效率、相關(guān)系數(shù)、相對誤差等。為了使模型能夠更好的模擬水沙過程，進(jìn)行調(diào)參數(shù)時優(yōu)先保證Nash效率、相關(guān)系數(shù)相對較高。

Nash-Sutcliffe效率：

Nash與Sutcliffe在1970年提出了模型效率系數(shù)（也稱確定性系數(shù)）來評價模型模擬結(jié)果的精度，它更直觀地體現(xiàn)了實測過程與模型模擬過程擬合程度的好壞，公式如下：

式中：Nash為Nash-Sutcliffe效率系數(shù)，其值越接近于1表示實測與模擬流量過程擬合得越好，模擬精度越高；為模擬值，為實測值，為實測平均值。

相關(guān)系數(shù)：

相關(guān)系數(shù)是對兩個變量之間關(guān)系的量度，考查兩個事物之間的關(guān)聯(lián)程度。相關(guān)系數(shù)的絕對值越大，相關(guān)性越強(qiáng)，相關(guān)系數(shù)越接近于1和－1，相關(guān)度越強(qiáng)，相關(guān)系數(shù)越接近于0，則相關(guān)度越弱。其計算公式如下：

式中：rxy為相關(guān)系數(shù)；n為系列的樣本數(shù)；X、Y分別代表實測系列和模擬系列的數(shù)值。通常情況下：|rxy|在0.8-1.0之間為極強(qiáng)相關(guān)，在0.6-0.8之間為強(qiáng)相關(guān)，在0.4-0.6之間為中等程度相關(guān)，在0.2-0.4之間為弱相關(guān)，在0-0.2之間為極弱相關(guān)或無相關(guān)。

相對誤差：

相對誤差是整個模擬期模擬值與實測值的差值與實測值的百分比，徑流量誤差絕對值越接近于零越好。

式中：Dv為模擬相對誤差（％）；F0為實測值均值；R為模擬均值。

到此為止，水沙耦合模擬模型已經(jīng)建立完成，如果參數(shù)經(jīng)過率定，就可以進(jìn)行水體保持計算了。

八、設(shè)置水土保持情景的步驟：用于使用水土保持措施數(shù)據(jù)，設(shè)置多個年代的水土保持措施情景。由于水沙耦合模擬模型是基于水沙物理過程建立起來的，其發(fā)生和發(fā)展過程通過輸入數(shù)據(jù)、計算公式進(jìn)行了確定了嚴(yán)格的計算關(guān)系。水土保持治理通過植樹、種草、封育、建設(shè)淤地壩等措施改變流域的水沙過程，對應(yīng)于模擬模型主要表現(xiàn)為土地利用、水利水保工程的變化。因此，通過改變模型用于計算的輸入數(shù)據(jù)和參數(shù)的變化，重新運(yùn)行模擬模型，即可以獲得不同水土保持治理條件下的流域水沙過程，通過比較和計算不同的模擬結(jié)果，定量評價不同水土保持治理取得的減水減沙效益。

九、生成不同水土保持情景下的水沙過程數(shù)據(jù)的步驟：用于將設(shè)計的水土保持情景數(shù)據(jù)運(yùn)行計算，得到多個水土保持措施的水沙過程數(shù)據(jù)。

十、分析計算不同水土保持情景的保水保土效益的步驟：用于對比不同水土保持情景下的水沙過程變化量進(jìn)行分析。

十一、比較并結(jié)束的步驟：用于對各個情景水土保持情景的水土保持?jǐn)?shù)據(jù)集所對應(yīng)的水沙過程時間段數(shù)據(jù)進(jìn)行比較，并輸出比較結(jié)果，結(jié)束計算過程。

匯集指定計算單元或流域出口的水沙數(shù)據(jù)作為一個水土保持?jǐn)?shù)據(jù)集所對應(yīng)的水沙過程時間段數(shù)據(jù)，并對多個水沙過程時間段數(shù)據(jù)進(jìn)行比較，并輸出比較結(jié)果，結(jié)束計算過程。由于坡面的水土保持隨著年代的變化也在不斷的變化，如：某一段時間（幾年）修建魚鱗坑，在過幾年在魚鱗坑的基礎(chǔ)上栽種樹木，在過幾年樹木被移去而建立梯田等，總之隨著時間的推移，坡面在不斷的變化，為此，可以按照不同的時間段所收集的坡面水土保持?jǐn)?shù)據(jù)，形成這一時段的水土保持?jǐn)?shù)據(jù)集，利用這些數(shù)據(jù)集分別計算各個時段的坡面水沙過程，并對這些不同時段的數(shù)據(jù)集所對應(yīng)的坡面水沙過程所形成的水沙過程時間段數(shù)據(jù)進(jìn)行比較，以獲取不同時段的水土保持措施所取得的效果。

水土保持情景實例：

利用構(gòu)建好的涇河流域分布式水沙耦合模型，以1956-2000年降雨情景等數(shù)據(jù)為基礎(chǔ)，分析了涇河流域人類活動（下墊面變化）影響下的流域水沙過程變化趨勢，以及梯田、壩地等水利水保工程措施減水減沙效益。

情景設(shè)置：

為了嘗試說明人類活動對流域產(chǎn)水產(chǎn)沙過程的影響，及水土保持工程措施水保效益等進(jìn)行水土保持治理過程中人們關(guān)心的問題，本實例共設(shè)定4種不同的人類活動情景：

情景1：歷史下墊面情景，即真實下墊面變化情景下模擬1956-2000年涇河流域水沙變化過程，用于還原真實下墊面變化情況下流域水沙變化情況。

情景2：1985年土地利用下墊面情景，即下墊面保持1985年土地利用和水土保持工程措施狀態(tài)不變條件下1956-2000年模擬涇河流域水沙變化過程。與情景1對比可以說明1985年水土保持措施對流域水沙關(guān)系的影響。

情景3：2000年土地利用下墊面情景，即下墊面保持2000年土地利用和水土保持工程措施狀態(tài)不變條件下1956-2000年模擬涇河流域水沙變化過程。與情景2對比可以說明，1985-2000年間人類活動減水減沙與增水增沙過程的對比關(guān)系。

情景4：2000年土地利用-水土保持措施下墊面情景，即下墊面保持2000年土地利用但不含水土保持工程措施狀態(tài)，模擬涇河流域1956-2000年水沙變化過程。與情景3對比可以說明不同降雨情景下梯田、壩地等工程類水土保持措施的水土保持效益。

使用本實施例建立的模型對人類活動對涇河流域水沙過程演變的影響分析：

研究采用代表時段方法，分別對1956-1959年、1960-1969年、1970-1979年、1980-1989年、1990-2000年、1956-1979年、1980-2000年等7個情景時段的河道水沙量變化進(jìn)行對比分析。

從表1不同下墊面情景河道徑流流量變化可以看出，涇河流域歷史下墊面情景（情景1）下1956-1959、1960-1969、1970-1979、1980-1989、1990-2000的流域平均年徑流量分別為17.91、24.07、19.96、17.35、8.20億m³。說明流域徑流量與降雨的變化趨勢一致，自60年代起流域徑流量開始明顯減少。根據(jù)徑流的變化可以以1980年為界，可以明顯的分為豐枯兩個時段，其中1956-1979的多年平均徑流量達(dá)21.33億m³，是1956-2000平均年徑流量的1.24倍；而1980-2000年的多年平均徑流量為12.56億m³，僅為1956-2000平均年徑流量的72.85%。

與歷史下墊面情景（情景1）相比，1985年土地利用下墊面（情景2）1956-1959、1960-1969、1970-1979、1980-1989、1990-2000年對應(yīng)情景時段的平均年流域徑流量分別變化2%、1%、1%、0%、-1%，變化幅度均在±2%以內(nèi)，說明情景2相對于情景1人類活動對流域的產(chǎn)匯流條件影響較小，但可以明顯看出1985年下墊面情況下，增大了1956-1979年的流域橫向水通量，說明該時期人類活動造成的增洪量大于減洪量。相對于情景1，情景2的1956-1979年多年平均徑流量增加0.27億m³，而1980-2000年年均徑流量減少0.07億m³，在1956-2000年時段總體上表現(xiàn)為徑流量增加；情景3在不同時段的徑流量均表現(xiàn)為增加趨勢，說明1956-2000年以來人類活動造成的徑流增量大于水土保持減水量。對比情景3和情景2可以發(fā)現(xiàn)，人類活動在扣除1985-2000年之間新增的水土保持措施減水效益的基礎(chǔ)上，在1956－1959、1960－1969、1970－1979、1980－1989、1990－2000時段等不同情況下仍然分別增加0.40、0.76、0.51、0.52、0.25 億m³。

2000年下墊面（情景3）及2000年下墊面減去工程水保措施（情景4）與情景1相比，在不同時段均表現(xiàn)出徑流增加趨勢。說明在不同時段，2000下墊面情景現(xiàn)的人類活動增水量大于水土保持措施的減水量。同時對比情景3和情景4可以看出，工程類水土保持措施減水效益明顯：1956-1959、1960-1969、1970-1979、1980-1989、1990-2000對應(yīng)時段分別減水0.76、1.09、1.3、0.74、0.98億m³，其中由于1956-1979時間段降雨較多，平均減水效益為1.13億m³，而在1980-2000年時間段為0.87億m³，總體上年均減水效益為1億m³。

從圖2中可以看出，涇河流域出口的徑流量與降雨量的年代變化增減變化趨勢一致，說明自然條件對涇河流域徑流的影響仍占主導(dǎo)地位。但不同情境下的響應(yīng)關(guān)系不同，其中情景4在1956-1959年到1960-1969年期間的徑流量增長率遠(yuǎn)大于降雨增長率，而從1970-1979年開始，在降雨減少速率較小的情況下徑流量迅速減少，說明情景4條件下對緩解氣候的變化能力最弱，其次為情景3和情景2。表明梯田、壩地等水土保持工程措施對調(diào)節(jié)流域降水-產(chǎn)流關(guān)系起到了重要作用。

表2可以看出，涇河流域歷史下墊面（情景1）下1956－1959、1960－1969、1970－1979、1980－1989、1990－2000年的平均年流域輸出沙量分別為2.22、2.53、1.94、2.07、0.62億t。相對于1956－1959年時段，1960－1969、1970－1979、1980－1989、1990－2000年時段的平均年輸出沙量分別變化14%、-13%、-7%、-72%。說明自60年代起流域輸出泥沙開始明顯減少，但與流域降雨、徑流量不同的是，在1980-1989年時段的侵蝕量較1970-1979年時段多出0.13億t，說明該時期的人類活動的增沙量大于減沙量。以1980年為界，流域在1956-1979年的多年平均輸沙量2.23億t，1980-2000年多年平均輸沙量1.31億t，分別是1956-2000年多年平均的1.24和0.73倍。

1）相比，1985年土地利用下墊面（情景2）1956－1959、1960－1969、1970－1979、1980－1989、1990－2000年對應(yīng)時段的平均年輸沙量分別變化2%、0.3%、1%、1%、-3%，說明在1956-1985年間人類活動的增沙量大于減沙量，在1990-2000年人類活動的增沙量也大于減沙量。這種情況在2000年下墊面（情景3）中體現(xiàn)更為直接，與歷史下墊面相比，分別變化5%、6%、5%、5%、3%，說明不同時期的人類活動增加的產(chǎn)沙量大于其減少的產(chǎn)沙量。對比情景3和情景2可以發(fā)現(xiàn)，在扣除1985-2000年之間新增的水土保持措減沙效益的基礎(chǔ)上，1956－1959、1960－1969、1970－1979、1980－1989、1990－2000年時段等不同情況下人類活動造成的流域增加產(chǎn)沙量分別為0.07、0.13、0.08、0.10、0.04億t。

對比情景3和情景4可以看出，工程類水土保持措施減沙效益明顯：1956－1959、1960－1969、1970－1979、1980－1989、1990－2000年對應(yīng)時段分別減沙0.23、0.32、0.28、0.16、0.12億t，其中1956-1979年時段平均減沙效益為0.29億t，而在1980-2000年時間段為0.14億t，總體上年均減沙效益為0.22億t。

從圖3，可以看出1956-1959、1960-1969和1970-1979年三個年代間的流域產(chǎn)沙量變化趨勢與降雨變化相同。而自1980-1989開始，水土保持措施等人類活動開始對流域減沙產(chǎn)生顯著的正面影響。其中情景4的減沙速率明顯較大，與情景2對比發(fā)現(xiàn)，自1980年以來的非工程水保措施對流域減沙的貢獻(xiàn)較大。

從圖4可以看出，單位徑流輸沙量中以情景4最大，而情景3最小。其中情景3在1956-1959、1960-1969、1970-1979年三個年代的單位徑流輸沙量顯著低于其他情景，說明1980年以來的水土保持措施使人類活動顯著改變了涇河流域的水沙關(guān)系，使流域加速侵蝕的狀況得到一定程度的緩解。

實施例二：

本實施例是實施例的改進(jìn)，是實施例一關(guān)于蒸發(fā)蒸騰計算的細(xì)化。本實施例所述的蒸發(fā)蒸騰計算：

計算單元內(nèi)（格柵或等高帶內(nèi)）的蒸發(fā)蒸騰包括來自植被濕潤葉面（植被截留水）、水域、土壤、城市地表面、城市建筑物等的蒸發(fā)，以及來自植被干燥葉面的蒸騰。計算單元的平均蒸發(fā)蒸騰量模型采用下式算出：

(1)

式中，F_W、F_U、F_SV、F_IR、F_NI分別為計算單元內(nèi)水域、不透水域、裸地-植被域、灌溉農(nóng)田及非灌溉農(nóng)田的面積率（%）；E_W、E_SV、E_U、E_IR、E_NI分別為計算單元內(nèi)水域、不透水域、裸地－植被域、灌溉農(nóng)田及非灌溉農(nóng)田的蒸發(fā)量或蒸發(fā)蒸騰量。

水域的蒸發(fā)量（E_w）由下述Penman公式算出：

(2)

式中，RN為凈放射量；G為傳入水中的熱通量；Δ為飽和水蒸氣壓對溫度的導(dǎo)數(shù)；δe為水蒸氣壓與飽和水蒸氣壓的差；r_a為蒸發(fā)表面的空氣動力學(xué)阻抗；ρ_a為空氣的密度；Ｃ_p為空氣的定壓比熱；λ為水的氣化潛熱，γ=Ｃ_p/λ。

裸地－植被域蒸發(fā)蒸騰量（E_SV）、灌溉農(nóng)田域（E_IR）、和非灌溉農(nóng)田域（E_NI）分別由以下公式計算:

(3)

E_IR=E_i1+E_i2+E_tr3+Es (4)

E_NI= E_i1+E_i2+E_tr4+Es (5)

式中，Ei為植被截留蒸發(fā)（來自濕潤葉面）；Etr為植被蒸騰（來自干燥葉面）；Es為裸地土壤蒸發(fā)。另外，下標(biāo)1表示高植被（森林、城市樹木），下標(biāo)2表示草，下標(biāo)3表示灌溉農(nóng)作物，下標(biāo)4表示非灌溉農(nóng)作物。

各類植被的截留蒸發(fā)（Ei）使用ISBA模型計算:

(6)

(7)

(8)

(9)

(10)

式中，Veg為植被覆蓋度；d為濕潤葉面積占總?cè)~面積的比例；E_p為可能蒸發(fā)量（由Penman方程式計算）；W_r為植被截留水量；P為降雨量；Rr為植被流出水量；W_rmax為最大植被截留水量；LAI為葉面積指數(shù)。

植被蒸騰由Penman－Monteith公式計算。

(11)

(12)

式中，RN為凈放射量；G為傳入植被體內(nèi)的熱通量；r_c為植物群落阻抗（canopy resistance）。蒸騰屬于土壤、植物、大氣連續(xù)體（SPAC：Soil－Plant－Atmosphere Continuum）水循環(huán)過程的一部分，受光合作用、大氣濕度、土壤水分等的制約。這些影響通過式（12）中的植物群落阻抗（r_c）來考慮。

植被蒸騰是通過根系吸水由土壤層供給。假定根系吸水率隨深度線性遞減、根系層上半部的吸水量占根系總吸水量的70%，則可得下式：

(13)

(14)

式中， Etr為蒸騰量；為根系層的厚度；z為離地表面的深度； S_r(z)為深度處的根系吸水強(qiáng)度；Etr (z)為從地表面到深度z處的根系吸水量。

根據(jù)以上公式，只要給出植物根系層厚，即可算出其從土壤層各層的吸水量（蒸騰量）。在本研究中，認(rèn)為草與農(nóng)作物等低植物的根系分布于土壤層的一、二層，而樹木等高植物的根系分布于土壤層的所有三層。結(jié)合土壤各層的水分移動模型，即可算出各層的蒸騰量。

裸地土壤蒸發(fā)由下述修正Penman公式計算：

(15)

(16)

式中，b為土壤濕潤函數(shù)或蒸發(fā)效率；q為表層（一層）土壤的體積含水率；q_fc為表層土壤的田間持水率；q_m為土壤單分子吸力（約1000~10000個大氣壓）對應(yīng)的土壤體積含水率。

不透水域的蒸發(fā)及地表徑流由下述方程式求解：

(17)

(18)

(19)

(20)

(21)

(22)

(23)

式中，P為降雨；H_u為洼地儲蓄；E_u為蒸發(fā)；R_u為表面徑流；H_umax為最大洼地儲蓄深；E_umax為潛在蒸發(fā)（由Penman公式計算）；c為城市建筑物在不透水域的面積率；下標(biāo)1表示城市建筑物，下標(biāo)2表示城市地表面。

實施例三：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于入滲計算的細(xì)化。本實施例所述的入滲計算采用Green-Ampt鉛直一維入滲模型模擬降雨入滲及超滲坡面徑流，和通用Green-Ampt模型進(jìn)行計算。

當(dāng)入滲濕潤鋒到達(dá)第m土壤層時入滲能力由下式計算：

(24)

式中，f為入滲能力；F為累積入滲量；k_m、A_m-1、B_m-1見后述。累積入滲量F的算出方法，視地表面有無積水而不同。

如果自入滲濕潤鋒進(jìn)入第m-1土壤層時起地表面就持續(xù)積水，那么累積入滲量由式（25）計算；如果前一時段t_n-1地表面無積水，而現(xiàn)時段t_n地表面開始積水，那么由式（26）計算。

(25)

(26)

(27)

(28)

(29)

(30)

(31)

式中，SW為入滲濕潤鋒處的毛管吸力；k為土壤層的導(dǎo)水系數(shù)；q_s為土壤層的含水率；q₀為土壤層的初期含水率；t為時刻；F_p為地表面積水時的累積入滲量；t_p為積水開始時刻；I_p為積水開始時的降雨強(qiáng)度；t_m-1為入滲濕潤鋒到達(dá)第m層與第m-1層交界面的時刻；L為入滲濕潤鋒離地表面的深度；L_i為第i層的厚度；Dq為q_s－q₀。

土壤水分吸力關(guān)系采用Havercamp公式：

(32)

式中，q為土壤體積含水率；q_s為飽和含水率；q_r為殘留含水率；j為吸引壓，cm水柱；a和b為常數(shù)。

土壤水分導(dǎo)水系數(shù)關(guān)系采用Mualem公式：

(33)

式中，K_s為土壤飽和導(dǎo)水系數(shù)，cm/s；k(q)為含水率q對應(yīng)的導(dǎo)水系數(shù)，cm/s；n為常數(shù)。

實施例四：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于地表徑流和壤中徑流計算的細(xì)化。本實施例所述的地表徑流和壤中徑流的計算。

地表徑流：

水域的地表徑流等于降雨減去降雨時的蒸發(fā)損失，不透水域的地表徑流按上述公式（20）及（22）計算，而裸地－植被域（透水域）的地表徑流則根據(jù)降雨強(qiáng)度是否超過土壤的入滲能力分以下兩種情況計算。

（1）霍頓坡面徑流（Hortonian overland flow）

當(dāng)降雨強(qiáng)度超過土壤的入滲能力時將產(chǎn)生這類地表徑流R1_ie，即超滲產(chǎn)流，由下式計算：

(34)

(35)

式中，P為降水量；H_SV為裸地－植生域的洼地儲蓄；H_SVMax為最大洼地儲蓄深；E_SV為蒸散發(fā)；f_SV為由通用Green-Ampt模型算出的土壤入滲能力。

（2）飽和坡面徑流（Saturation overland flow）

對于河道兩岸及低洼的地方，由于地形的作用，土壤水及淺層地下水逐漸匯集到這些地方，土壤飽和或接近飽和狀態(tài)后遇到降雨便形成飽和坡面徑流（蓄滿產(chǎn)流）。此時，Green-Ampt模型已無能為力，需根據(jù)非飽和土壤水運(yùn)動的Richards方程來求解。為減輕計算負(fù)擔(dān)，地表洼地儲留層按連續(xù)方程、表層土壤分成3層按照Richards方程（積分形式）進(jìn)行計算：

地表洼地儲留層

(36)

(37)

土壤表層

(38)

土壤中層

(39)

土壤底層

(40)

( j：1、3) (41)

(42)

(43)

( j：1、2) (44)

( j：1、2) (45)

式中，H_s為洼地儲蓄；H_smax為最大洼地儲蓄；Veg₁、Veg₂為裸地－植生域的高植生和低植生的面積率；Rr₁、Rr₂為從高植生和低植生的葉面流向地表面的水量；Q為重力排水；QD_j,j+1為吸引圧引起的j 層與j+1 層土壤間的水分?jǐn)U散；E₀為洼地儲蓄蒸發(fā)；E_s為表層土壌蒸發(fā)；E_tr為植被蒸散（第一下標(biāo)中的1表示高植生、2表示低植生；第一下標(biāo)表示土壤層號）；R₂為壤中徑流；k(q)為體積含水率q對應(yīng)的土壤導(dǎo)水系數(shù)；y(q)為體積含水率q對應(yīng)的土壤吸引力；d為土壤層厚度；W為土壤的蓄水量，w=qd）；W₁₀為表層土壤的初期蓄水量。另外、下標(biāo)0、1、2、3分別表示洼地儲蓄層、表層土壤層、第2土壤層和第3土壤層。

壤中徑流：

在山地丘陵等地形起伏地區(qū)，同時考慮坡向壤中徑流及土壤滲透系數(shù)的各向變異性。從不飽和土壤層流入河道的壤中徑流由下式計算：

(46)

式中，k(θ)為體積含水率q對應(yīng)的沿山坡方向的土壤導(dǎo)水系數(shù)；slope為地表面坡度；L為計算單元內(nèi)的河道長度；d為不飽和土壤層的厚度。

實施例五：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于坡面匯流計算的細(xì)化。本實施例所述的坡面匯流計算為：

坡面匯流。模型則在河網(wǎng)水系生成和流域劃分的基礎(chǔ)上，根據(jù)網(wǎng)格單元DEM及土地利用等基本信息準(zhǔn)備各子流域等高帶的屬性表（包括面積、長度、寬度、平均高程、坡度和曼寧糙率等），采用一維運(yùn)動波法從上游等高帶至下游等高帶計算坡面匯流，并將下游等高帶的坡面匯流輸入給所在子流域內(nèi)的河道。

運(yùn)動波方程：

（連續(xù)方程） (50)

（運(yùn)動方程） (51)

（Manning公式） (52)

式中，A為流水?dāng)嗝婷娣e；Q為流量；q_L為網(wǎng)格單元或河道的單寬流入量（包含網(wǎng)格內(nèi)的有效降雨量、來自周邊網(wǎng)格及支流的水量）；n為Manning糙率系數(shù)；R為水力半徑；S₀為網(wǎng)格單元地表面坡降或河道的縱向坡降；S_f為摩擦坡降。

動力波方程（Saint Venant 方程）：

（連續(xù)方程） （53)

（運(yùn)動方程） (54)

（Manning公式） (55)

式中，A為流水?dāng)嗝婷娣e；Q為斷面流量；q_L為網(wǎng)格單元或河道的單寬流入量（包含網(wǎng)格內(nèi)的有效降雨量、來自周邊網(wǎng)格及支流的水量）；n為Manning糙率系數(shù)；R為水力半徑；S₀為網(wǎng)格單元地表面坡降或河道的縱向坡降；S_f為摩擦坡降；V為斷面流速；V_x為單寬流入量的流速在x方向的分量。

實施例六：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于積雪融雪過程計算的細(xì)化。本實施例所述的積雪融雪過程計算：

（56）

（57）

式中，SM為融雪量，mm/day；M_f為融化系數(shù)或稱“度日因子”，mm/°C/day；T_a為氣溫指標(biāo)，°C；T₀為融化臨界溫度，°C；S為積雪水當(dāng)量，mm；SW為降雪水當(dāng)量，mm；E_snow為積雪升華量，mm。

“度日因子”既隨海拔高度和季節(jié)變化，又隨下墊面條件變化，常做為模型調(diào)試參數(shù)對待。一般情況下在1 mm/°C/day和7 mm/°C/day之間，且裸地高于草地，草地高于森林。氣溫指標(biāo)通常取為日平均氣溫。融化臨界溫度通常在-3°C和0°C之間。另外，為將降雪與降雨分離，還需要雨雪臨界溫度參數(shù)（通常在0°C和3°C之間）。

實施例七：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于雨滴濺蝕模擬和計算的細(xì)化。本實施例所述的雨滴濺蝕模擬和計算：

（58）

式中，D₁為雨滴擊濺侵蝕量，g/m²；E_drop為雨滴動能，J/m²；I為雨強(qiáng)，mm/min；J₁為地表坡度，°；k₁，α₁，β₁為經(jīng)驗參數(shù)。

其中雨動能E的計算：

（59）

式中，E_unit為單位降雨動能，J/(m²·mm)；k₁＇,α₁＇為經(jīng)驗參數(shù)；

當(dāng)水深大于雨滴直徑3倍以上時，取水深大于0.6cm時，雨滴濺蝕作用消失；

由雨滴濺蝕增加的土壤侵蝕輸沙能力計算公式為：

（60）

式中，q_s1為單寬流量輸沙能力，kg/m²；k₂，α₂，β₂為經(jīng)驗參數(shù)。

實施例八：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于薄層水流侵蝕模擬的細(xì)化。本實施例所述的薄層水流侵蝕模擬：

（61）

（62）

（63）

式中，Dc為水流剝離土壤速率；k₃為土壤的可蝕性參數(shù)；τ_f為水流對土壤顆粒的剪切應(yīng)力；τ_c為土壤的臨界抗剪切應(yīng)力；Dr為細(xì)溝水流剝蝕率；q為單寬流量；c為水流泥沙含量； Tc為水流的挾沙能力；k₄，α₄為經(jīng)驗常數(shù)。

實施例九：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于股流侵蝕過程模擬的細(xì)化。本實施例所述的股流侵蝕過程模擬，股流的挾沙能力T_SE計算公式為：

（64）

式中，k₅為淺溝水流攜沙能力系數(shù)；m為側(cè)向匯流影響常數(shù)，ω_u為單位水流功率。

（65）

式中，S_E為股流侵蝕量；Q_E為流量；D_r為上游來沙量。

實施例十：

本實施例是上述實施例的改進(jìn)，是上述實施例關(guān)于重力侵蝕過程模擬的細(xì)化。本實施例所述的重力侵蝕過程模擬，重力侵蝕量計算公式為：

坡面水沙過程中的重力侵蝕，其影響因素相對單一。對于坡面坡度大于黃土休止角的土體，其是否發(fā)生重力侵蝕取決于土體的下滑力與土體抗剪強(qiáng)度之間的大小關(guān)系。當(dāng)重力侵蝕發(fā)生。如圖2所示，土體下滑力可用下式計算：

（72）

式中， Δh為發(fā)生重力侵蝕土體深度，細(xì)溝溝坡取0.1m，淺溝溝坡取0.5m，切溝溝坡取10m；θ為黃土天然休止角；γ_s為土壤濕容重。

重力侵蝕量計算公式為：

（73）

式中，為發(fā)生重力侵蝕的溝長系數(shù)，為溝道總長度，其它參數(shù)意義同上。

土壤的臨界抗剪強(qiáng)度公式為：

，n=78，r=0.61 （74）

由于坡面侵蝕過程中重力侵蝕發(fā)生的尺度相對較小，影響因素相對單一。采用土體的下滑力與土體的抗剪強(qiáng)度大小作為判定條件。其中土體的下滑力由式（72）計算，土體的抗剪強(qiáng)度由式（74）計算，重力侵蝕量由式（73）計算。

最后應(yīng)說明的是，以上僅用以說明本發(fā)明的技術(shù)方案而非限制，盡管參照較佳布置方案對本發(fā)明進(jìn)行了詳細(xì)說明，本領(lǐng)域的普通技術(shù)人員應(yīng)當(dāng)理解，可以對本發(fā)明的技術(shù)方案（比如公式的運(yùn)用、步驟的先后順序等）進(jìn)行修改或者等同替換，而不脫離本發(fā)明技術(shù)方案的精神和范圍。