2新安江流域水文模型

2024年1月12日發(作者：三國演義歷史背景)

2新安江流域水文模型60年代初，河海大學(原華東水利學院)水文系趙人俊等開始研究蓄滿產流模型，配合一定的匯流計算，將模型應用于水文預報和水文設計。1973年，他們在對新安江水庫做入庫流量預報的工作中，把他們的經驗歸納成一個完整的降雨徑流流域模型——新安江模型。模型可用于濕潤地區和半濕潤地區的濕潤季節徑流模擬和計算。最初的新安江模型為兩水源模型，只能模擬地表徑流和地下徑流。80年代初期，模型研制者將薩克拉門托模型與水箱模型中，用線性水庫函數劃分水源的概念引入新安江模型，提出了三水源新安江模型，模型可以模擬地面徑流、壤中流、地下徑流。1984至1986年，又提出了四水源新安江模型，可以模擬地面徑流、壤中流、快速地下徑流和慢速地下徑流。三水源新安江模型一般應用效果較好，但模擬地下水豐富地區的日徑流過程精度不夠理想。在新安江三模型中增加慢速地下水結構就成為四水源新安江模型。 當流域面積較小時，新安江模型采用集總模型，當面積較大時，采用分塊模型。分塊模型把流域分成許多塊單元流域，對每個單元流域做產、匯計算，得到單元流域的出口流量過程。再進行出口以下的河道洪水演算，求得流域出口的流量過程。把每個單元流域的出流過程相加，就求得了流域出口的總出流過程。劃分單元流域的主要目的是處理降雨分布的不均勻性，因此單元流域應當大小適當，使得每塊面積上的降雨分布比較均勻。并有一定數目的雨量站。其次盡可能使單元流域與自然流域相一致，以便于分析與處理問題，并便于利用已有的小流域水文資料。如果流域內有大中型水庫，則水庫以上的集水面積即應作為一個單元流域。因為各單元流域的產匯、流計算方法基本相同，以下只討論一個單元流域的情況。新安江模型包括4個計算環節：蒸散發計算；流域產流計算；徑流

劃分；匯流計算。4個計算環節分別概化了流域降雨徑流的主要產、匯流物理過程。2.1流域蒸散發計算各種水源的蒸散發計算模型均可采用兩層蒸發模型或兩層蒸發模型，一般根據實際情況選用。原則是在模擬徑流精度相同的情況下，盡量采用參數少的兩層蒸散發模型。蒸散發模型不考慮面上分布的不均勻性，但可考慮土濕垂向分布的不均勻性。2.2.1兩層蒸散發模型兩層蒸散發模型將土層分為上、下兩層，各層蓄水容量分別為、（）。流域土層下滲蓄水和蒸散發計算過程按下述原則進行：降雨先補充上層，上層蓄滿后再補充下層；蒸發時先蒸發上層，上層蓄水量蒸發殆盡后再蒸發下層。兩層蒸散發模型的計算式如下：⑴：， ，⑵ ： ，，式中：——是流域降雨量，；——是流域蒸散發能力，；——是上層土壤蓄水量，；——是下層土壤蓄水量，；——是下層土壤蓄水容量，；——上層蒸發量，?！聦诱舭l量，。2.2.2三層蒸散發模型三層蒸散發模型將土層分為上層、下層、深層三個土層，各層蓄

水容量分別為、、（）。流域土層下滲蓄水和蒸散發的計算過程與兩層蒸發模型類似，按下述原則進行：降雨首先補充上層，上層蓄滿后再補充下層，下層蓄滿后再補充深層；蒸發時先蒸發上層，上層蓄水量蒸發殆盡后再蒸發下層，下層蓄水量蒸發殆盡后再蒸發深層。由土壤蒸發規律知，在久旱時，流域下層土壤含水量趨于干涸，深層土壤水分將以薄膜水、氣態水的形式向上輸送，使蒸發維持在一個較小的穩定數值。因此，考慮來自了深層的水分蒸發。為了使流域的蒸散發不小于某個最小的穩定蒸散發量，模型中考慮了一個深層蒸散系數，并在結構上進行了控制，使計算過程中的最小蒸散發量不會小于。三層蒸散發模型的計算式如下：⑴當：， ， ，；⑵ ，時：，，，；⑶，時：，，，；⑷，時：，，，。式中：——是深層蒸散發系數。2.2新安江兩水源模型模擬流域出口斷面的徑流過程時，如果水源劃分只考慮直接徑流和地下徑流，而不考慮壤中流，這樣的新安江流域模型稱為新安江兩水源模型。2.2.1模型的結構新安江兩水源模型的4個計算環節中（蒸散發計算、流域產流計算、徑流劃分、匯流計算），蒸發計算可采用兩層蒸發模型，也可采用三層蒸發模型。產流計算采用蓄滿產流模型。徑流劃分計算利用穩定下滲率

將徑流劃分為地面徑流（或稱直接徑流）和地下徑流兩種水源。地面徑流采用單位線匯流，地下徑流采用一次線性水庫匯流，即地下匯流計算的坡面匯流和河網匯流用一個線性水庫的出流過程概化。模型把流域面積劃分為透水面積和不透水面積兩部分，不透水面積上的降水在滿足蒸發后將直接轉化為地面徑流。透水面積上將發生下滲，下滲的水量一部分存儲于土層（降雨過后耗于蒸發）。滿足了流域土壤蓄水容量后的下滲水量才能轉化為徑流。不透水面積用參數表示，它是用流域內不透水面積占全流域面積的百分比表示的。新安江模型的輸出是流域出流過程和流域蒸散發過程，輸入則為時段降雨量、蒸發皿觀測蒸發量。新安江兩水源模型結構如下圖所示：圖2.1-1 新安江兩水源模型結構示意圖圖中方框內標注為狀態變量，方框外標注的是模型參數。

2.2.2模型的產流計算新安江模型產流部分的計算是蓄滿產流模式， 蓄滿產流指在流域包氣帶土壤濕度達到田間持水量以前不產流，所有的降雨都被土壤吸收；而在土壤濕度達到田間持水量之后，所有的降雨（除去同期的蒸散發）都產流。在產流后，流域包氣帶土壤的下滲能力為穩定下滲率，下滲的水分成為地下徑流和壤中流，超蓄的部分成為地面徑流?？紤]到流域內各點的蓄水容量并不相同，實際產流常常是在部分面積上產圖2.2-1流域蓄水容量曲線示意圖 流，新安江模型引入一條流域蓄水容量曲線來刻劃流域內各點蓄水容量的不均勻性，把流域內各點的蓄水容量概化成如圖2.3所示的一條拋物線（也可概化成其它函數形式），其方程為： （2.2-1）式中：為產流面積；——是點蓄水容量, ；——是蓄水容量小于等于蓄水量的全部點的面積,；——是單元流域面積, ；——是蓄水容量曲線指數。——是流域內點最大蓄水容量,；其與流域平均蓄水容量的關系為: （2.2-2）引入蓄水容量曲線后，在降雨過程中，只有在滿足了包氣帶蓄水容量的面積上才可能產生徑流，其余面積上不產生徑流?？鄢舭l后的有效降雨量中，未轉化為凈雨產流量的水分通過下滲進入流域包氣帶土層，補充土層蓄水量。引入流域蓄水容量曲線后的產流計算式為：

當時為部分面積產流

（2.2-3）當時為全流域面積產流

（2.2-4） 式中：

——為有效降雨量， ；——為與流域土壤初始蓄水量相應的前期影響雨量，由下式計算。 （2.2-5）2.2.3模型的水源劃分及匯流計算1.新安江兩水源模型的水源劃分新安江兩水源模型的水源劃分只考慮直接徑流和地下徑流，在計算出流域總產流量后，進一步利用穩定下滲率將總徑流量劃分為地面徑流（或稱直接徑流）量和地下徑流量兩種水源。兩水源劃分的計算式如下： (2.2-6)

(2.2-7)2.新安江兩水源模型的匯流計算直接徑流采用一條單位線匯流直接匯流到流域出口斷面，即：

(2.2-8)式中：——第時刻的直接徑流流量，；——第時段的直接徑流量，；——單位線第時段末的流量，。地下徑流采用一個線性水庫匯流到流域出口斷面，即： (2.2-9)

式中：——單位轉換系數，，是流域面積，。流域出口斷面時刻的流量為： (2.2-10)2.2.4模型的參數及調試 采用三層蒸發計算的新安江兩水源模型共有9個參數，一條單位線。 ——流域蒸發折算系數，是流域蒸散發能力與蒸發皿蒸發量之比； ——深層蒸散發系數； ——不透水面積占全流域面積的百分比，％； ——流域平均蓄水容量(指張力水)，； ——流域平均上層蓄水容量, ； ——流域平均下層蓄水容量, ；——蓄水容量曲線指數；——穩定下滲率；——地下水消退系數；——單位線。當模型初定后就可上機調試，在調試的過程中，應注意各參數的徑流響應特征，以及參數之間相依性的影響，以便調試能有目的地順利進行。以下簡介的調試技術可供參考。流域蒸發折算系數的調試，對一個具體流域來說，此參數完全靠優選。調試此參數主要考察年徑流模擬誤差，當年徑流模擬誤差達最小時，此值最優，但個別參數的最優并不表示模型總體最優。根據現有經驗，年徑流模擬誤差控制在5％～8％左右就可以了。流域平均蓄水容量的調試，一般根據實測的降雨徑流資料分析或經驗選取后，適當微調即可。經驗選定后，此值調試范圍不大，比較容易確定。通常在南方濕潤地區可選80～150左右，江淮一帶約在110左右,

燕山東北東部地區，選150左右。流域平均上層蓄水容量，一般取10～20，流域平均下層蓄水容量一般取60～90。蓄水容量曲線指數的調試，表示流域蓄水容量分布的不均勻性，當全流域蓄滿后，這個參數就不起作用了，因此，應該選取流域沒有達到蓄滿的那些洪水點據作為調試的依據。此值一般在0.2～0.5之間取值。注意值與值之間有相依性，二者對徑流模擬的結果相互有影響。穩定下滲率的調試，對水源劃分起決定性作用，但至今對的研究還不充分，目前可對作一些簡單處理，如果取為一個常數使模擬效果不好時，可考慮作為變動參數?？疾焓欠窈线m主要觀察地下徑流的模擬精度。地下水消退系數的調試，主要觀察洪水退水段的擬合精度。不透水面積參數對模型擬合精度影響很小，一般取0～0.05。2.3新安江三水源模型模擬流域出口斷面的徑流過程時，如果水源劃分考慮直接徑流、壤中流和地下徑流三種水源，這樣的新安江流域模型稱為新安江三水源模型。2.3.1新安江三水源模型結構新安江三水源的4個計算環節中（蒸散發計算、流域產流計算、徑流劃分、匯流計算）。流域的蒸散發計算和產流計算與新安江二水源模型相同，水源劃分則完全不同，因此匯流計算也不相同。新安江三水源模型結構示意圖如下圖所示。圖2.3-1 新安江三水源模型結構示意圖2.3.2新安江三水源模型的水源劃分

1.不考慮自由水蓄水容量分布不均勻的劃分方法求得的產流量包括地面徑流 、壤中流和地下徑流 三部分。新安江三水源模型用一個自由水蓄水庫解決水源劃分，自由水蓄水庫有兩個出流孔，底孔為地下徑流出流孔，邊孔為壤中流出流孔。新安江模型考慮了產流面積（）的變化，自由水蓄水庫實際只發生在產流面積上，其底寬為產流面積，顯然它是隨時間變化的。產流量進入水庫即在產流面積上產生的徑流深，也就是自由水蓄水庫所增加的蓄水深，當自由水蓄水深超過其最大值時，超過部分成為地面徑流。壤中流和地下徑流按線性水庫出流，其出流系數分別為和。底孔出流量和邊孔出流量分別進入各自的水庫，并按線性水庫的退水規律流出，分別成為地下水出流和壤中流出流。模型認為蒸散發在張力水中消耗，自由水蓄水庫的水量全部為徑流。水源劃分示意圖見圖2.3-2。水源劃分計算式為：圖2.3-2 不考慮自由水蓄水容量分布不均勻的劃分示意圖地面徑流計算式：當 (2.3-1)當 (2.3-2)壤中流計算式：

(2.3-3)地下徑流計算式：

(2.3-4)注意每步計算時應將自由水蓄水庫的蓄量折算為當前時段產流面積上的深度。2.考慮自由水蓄水容量分布不均勻的劃分方法流域自由水蓄水容量實際上也是不均勻的，作為圖2.3-1表示的水源劃分方法的改進，人們引入流域自由水蓄水容量分布曲線來刻化自由水蓄水量的不均勻性，其線型與流域蓄水容量曲線類似，曲線見圖2.3-3，曲線方程如下式。 （2.3-5）

式中：為產流面積；——是點自由水蓄水容量；——是自由水蓄水容量小于等于蓄水量的全部點的面積；——是單元流域面積；——是自由水蓄水容量曲線指數?！橇饔騼赛c最大自由水蓄水容量；與流域平均自由水蓄水容量的關系為： （2.3-6）由于認為在產流面積上才有自由水，因此，產流面積上的點最大自由水蓄水容量是隨產流面積變化而變化的，見圖2.3-3。對當前時段，已知，則將其代人（2.3-5）得 （2.3-7）1.0圖2.3-3 考慮自由水蓄水容量分布不均勻的劃分示意圖從（2.3-7）解出得： （2.3-8）同樣，在產流面積上的自由水平均蓄水容量與該面積上的點最大自由水蓄水容量的關系為： （2.3-9）在劃分水源時，產流面積上只有在滿足了蓄水容量的水庫面積上才可能產生地面徑流；其余面積上不產生地面徑流。進入水庫的水分扣

除地面徑流量后，剩余的補充水庫蓄水量S，這時的水源劃分計算式為：⑴地面徑流出流當時為部分面積產流

（2.3-10）當時為全流域面積產流

（2.3-11）式中：——為計算時段初的自由水蓄水庫的蓄水量（用產流面積上的平均深度表示）；——為計算相應于的自由水蓄水庫前期影響入流量，由下式計算； （2.3-12）進入水庫的水分扣除地面徑流量后，剩余的補充水庫蓄水量： （2.3-13）自由水蓄水庫的時段入流量為時段產流量，但由（2.3-3）和（2.2-4）計算的產流量是用流域平均水深表示的水量，應將折算為產流面積上表示的水深?？梢宰C明折算為產流面積上表示的水深值恰為流域有效降雨值，證明如下：為折算為產流面積上表示的水深，根據降雨～徑流關系的性質，水文學原理中已經證明相對產流面積與產流量之間有關系：

（2.3-14）從此式可看出⑵壤中流出流：

（2.3-15）⑶地下徑流出流： （2.3-16）由于,都是產流面積上的深度表示的水量，計算所得、、應折算回

全流域深度表示的水量，即：、、在當前時段計算結束后應乘。2.3.3流域匯流計算匯流計算可以采用兩種方法：1.傳統流域匯流這種匯流方法是：地面徑流采用單位線匯流，壤中流和地下徑流采用線性水庫匯流至出口斷面后與單位線匯流結果疊加，最后得到流域出口斷面的模擬出流。即圖2.3-2中下面的兩個水庫的出流就是壤中流和地下徑流在流域出口斷面的出流流量，該兩個水庫分別稱為壤中流調蓄水庫和地下徑流調蓄水庫。⑴壤中流匯流計算 （2.3-17）⑵地下徑流匯流計算 （2.3-18）式中： ——時段壤中流出流流量，；——壤中流消退系數，可假定為；——時段地下徑流出流流量，；——地下徑流消退系數，可假定為；——是徑流深轉換為流量時的單位轉換系數，.上式中假定了調蓄水庫與自由水水庫的出流系數相同。在實踐中，為延緩壤中流和地下徑流出流，有人把下面的兩個水庫的出流看作是坡地調蓄水庫，而在它們后面再分別串聯一個水庫反映壤中流和地下徑流受河槽調蓄的作用。2.河網匯流單位線匯流河網匯流單位線也稱時變河網匯流單位線，這種匯流方法忽略地面徑流的坡面匯流階段，直接將地面徑流產流量轉換為流量作為河網地面徑流入流（看作集中在計算時段末入流）。地下徑流和壤中流匯流均采

用線性水庫作坡地調蓄后得到相應河網入流、，即：把圖2.3-2或圖2.3-3中的轉換為流量作為河網地面徑流入流，下面的兩個水庫的出流看作壤中流和地下徑流經坡地調蓄后注入河網的入流量。三種水源的河網入流流量相加就是所謂的河網總入流。河網總入流經流域河網調蓄后得到流域出口斷面的流量過程。河網匯流單位線的計算式如下：

（2.3-19）式中： ——時段河網的總入流，，；——河網消退系數；——河網匯流時間，；——時刻的河網出流（即時刻的流域出口斷面流量），；——時刻的河網總入流經河網調蓄后，延遲時間的出流，。是隨時間而變化的時變參數，它取決于時刻的河網的總入流量和流域的河槽特性，用下式計算： （2.3-20）式中： ——反映流域河槽特性的系數。2.3.4線性水庫不同單位時段的出流系數轉換 從(2.3-17)和(2.3-18)，知道新安江模型的自由水水庫出流是按線性水庫出流處理的。在客觀世界中，流域的降雨、蒸散發、徑流等水文物理量都是連續量，流域內的蓄水量也是連續量，觀測和計算時都只能處理為一個一個的時段平均量；從數學觀點來看，這實際上是將連續量處理為離散量，將表達物理過程的連續方程處理為表達離散量的差分方程，這樣就必然帶來差分誤差，特別是水文模型計算通常采用的時段較大，如果在模型的算法設計上不注意，差分誤差可能較大。為了減小計算的差分誤差，本可以縮小計算時段，但是，在水文上常常受到限制，這是由于原始水文觀測數據時段本身較大，如果在所有的計算環節上都縮小計算時段并不能提高成果的精度。盡管如此，在模型的算法上，至

少應要求算法不會帶來更大的差分誤差。仔細考察三水源新安江流域水文模型計算的各個環節，就會發現(2.3-17)和(2.3-18)的自由水水庫出流計算可能會出現較大的差分誤差，為此，有的模型程序設計者在算法上作了一定處理，以設法減小差分誤差?；镜南敕ㄊ窃谶@個計算環節上減小計算步長，即將模型總體的計算時段在這個環節上進一步細分；為此，需要將水庫長時段的出流系數轉換為短時段的出流系數。設水庫在一個計算時段內的入流量為，為減小差分誤差，將等分為份入流，這實際上是將原計算時段細分為個更小的計算時段，每個時段的水庫入流為?？梢宰C明水庫計算時段為的出流系數和與計算時段為的出流系數和之間存如下關系： ； （2.3-21）證明如下：1.只有一個出流孔的線性水庫不同時段出流系數的轉換先對只有一個出流孔的線性水庫進行證明。設線性水庫的出流方程和退水時的水量平衡方程如下：

（2.3-22）

（2.3-23）式中： ——為線性水庫的瞬時出流系數；——為線性水庫時刻的蓄量；——為線性水庫時刻的出流流量。合并（2.3-22）和（2.3-23）得到: （2.3-24）在上積分上式得 （2.3-25）此式稱為線性水庫退水方程。此式的另一個表達形式是：

（2.3-26）這是（2-27）在區間上積分的結果。令，即，取為一個單位計算時段，則上式記為： （2.3-27）稱為線性水庫消退系數（又稱退水系數）。在區間上積分，得到水庫在時段內流出的水量。

（2.3-28）令，即，取為一個單位計算時段，則上式記為 （2.3-29）（2.3-29）中的就是線性水庫單位時段為的出流系數。比較（2.3-29）和（2.3-22）可知和的意義是完全不同的，總是與一個單位時段相聯系的，單位時段不同，其值也不同。令：

記單位時段為的水庫消退系數為。由（2.3-27）的退水方程，可推導得 （2.3-30）此式中的1代表一個時段，表示個時段。顯然。對單位時段的退水方程，有； 或 （2.3-31）這里，1代表一個時段。比較（2.3-30）最后一個式子與（2.3-31）式可知，線性水庫單位計算時段為的出流系數與單位計算時段為的出流系數之間有如下關系：

（2.3-32）已知對單位計算時段，有；記單位計算時段的出流系數為，則可推得 （2.3-33）此式就是所需要得結果。2.有二個出流孔的線性水庫不同時段出流系數的轉換記單位計算時段為的出流系數分別為和，單位計算時段為的出流系數分別為和，不同計算時段的兩孔總出流系數分別為和，代人（2.3-

33）式得到： （2.3-34）要求解和還須增加一個條件，可假定不同計算時段的兩孔出流系數分配比例相同，即成立： （2.3-35）利用（2.3-34）和（2.3-35）消去得到： （2.3-36）利用（2.3-34）和（2.3-35）消去得到： （2.3-37）利用2.3-36）和（2.3-37）就可進行線性水庫不同時段的出流系數轉換。2.3.5新安江三水源模型的參數調試采用河網匯流單位線的新安江三水源模型共有14個參數，根據參數的功能可把參數分為4類： ⑴與蒸散發計算有關的參數：蒸散發折算系數，上層土壤蓄水容量下層土壤蓄水容量，深層蒸散發系數。⑵與產流計算有關的參數：流域蓄水容量，流域蓄水容量曲線指數，流域不透水面積比值。⑶與水源劃分有關的參數：流域自由水蓄水容量，流域自由水蓄水容量曲線指數，自由水蓄水庫壤中流出流系數與自由水蓄水庫地下徑流出流系數。⑷與匯流有關的參數：壤中流消退系數和地下徑流消退系數；河網匯流單位線參數和。壤中流匯流調蓄水庫和地下徑流匯流調蓄水庫與自由水水庫是不同的水庫，但如果假定二調蓄水庫的出流系數與自由水水庫對應的徑流的出流系數相同，即假定，，則參數可減為12個。4類參數相對獨立，同類參數有一定相依性。1.參數初定根據流域自然地理資料、水文氣象資料和前人以積累的經驗初定參

數。⑴蒸散發折算系數可采用汛期連續大雨，流域蓄水量處于蓄滿狀態的累計降雨、累計徑流、累計蒸發觀測資料，用下式估計: （2.3-38）⑵流域蓄水容量、上層土壤蓄水容量、下層土壤蓄水容量值一般根據實測的降雨徑流資料分析或經驗選取后，適當微調即可。經驗選定后，此值調試范圍不大，比較容易確定。通常在南方濕潤地區可選80～150左右，北方半濕潤地區可達170，江淮一帶約在110左右,燕山東北東部地區，選150左右。 中包含有植物截留，少林地5～10，多林地10～20左右。在60～90范圍。 ⑶流域蓄水容量曲線指數此值一般在0.2～0.5之間取值，通常流域越小，該值越小，特小流域可能小于0.2。1000平方公里以上流域一般大于0.3。 ⑷不透水面積比值值對模型擬合精度影響很小，一般取0～0.05。⑸深層蒸散發系數 南方多林地區左右 ；北方半干旱地區左右。 ⑹流域自由水蓄水容量 一般由經驗初定，土層淺薄的山區，土厚多林地區，甚至更大。 ⑺流域自由水蓄水容量曲線指數 根據經驗，常在1～1.5左右。⑻自由水蓄水庫壤中流出流系數與自由水蓄水庫地下徑流出流系數一般來說，約為0.7左右，二者分配比與具體流域有關，只能調試決定。⑼壤中流消退系數和地下徑流消退系數：一般趨近0.9；約為0.98～0.998，具體取值與具體流域有關，可

進一步調試決定。⑽河網匯流單位線參數和在0.02左右，與具體流域的河網特性有關，可進一步調試決定。與具體流域的河網特性有關，為計算方便，通常取計算單位時段的整倍數。2.參數率定首先率定第一類參數，主要觀察率定期的總水量是否與實測水量相近。其次率定第二類參數，主要觀察次洪水的水量是否與實測水量相近。第三類參數的率定，主要觀察次洪水的退水段過程的水量是否與實測水量相近。第四類參數的率定，主要觀察次洪水的過程線是否擬合良好。參數率定過程中要反復對參數進行多遍調試，直至滿足精度要求。2.3.6新安江三水源模型的VB程序介紹新安江三水源，兩層蒸發模型（時變單位線）VB源程序Function XingAnJang4() As Integer Dim Fno As Integer, I As Integer，Nian(2) As Integer，AX(20)As Single, Qo As Single Dim Wm As Single, Wum As Single, Wlm As Single, k As Single,Sm As Single, b As Single Dim Bx As Single, FR0 As Single, FAs Single, U As Single, X As Single, XX As Single Dim Kss As Single, Kg As Single, CRss As Single, CRg AsSingle, Immp As SingleDim P As Single, E As Single, Pe As Single, Eu As Single, El

As Single, EE As Single

Dim R As Single, Rs As Single, Rss As Single, Rg As Single, QAs Single, Qz As Single,Dim W As Single, Wu As Single, Wl As Single, S As Single, WmmAs Single, Smm As Single

Dim QRss0 As Single, QRg0 As Single, Cr As Single, Cs AsSingle, Cx As SingleDim Qa() As Single Qo = 5 'Qo是基流。 Open "e:test流域模型" For Input As #1 Open "e:test流域模型" For Input As

#5 Open "e:test流域模型" For Output As

#6 Fno = 1ReadRecord Fno:TurnRecord Nian, AX, 11 + 1, 1: Clo (Fno)Wm = AX(1): Wum = AX(2): k = AX(3): b = AX(4)Sm = AX(5): Bx = AX(6): Kss = AX(7) Kg = AX(8)CRss = AX(9): CRg = AX(10): Immp = AX(11) W = 85: Wu = 1: S = 0: FR0 = 0: Dt = 3: F = 2928: U = F /(3.6 * Dt): QRss0 = 0 QRg0 = 10.1 Wmm = (1 + b) * Wm: Smm = (1 + Bx) * Sm: Wlm = Wm - WumWl = W - Wu: X = 1 / (1 + b): XX = 1 / (1 + Bx) Cr = 0.025: Cx = 6

Fno = 5: n = 0: Do Until EOF(Fno): Line Input #Fno, TextLien:n = n + 1: Loop: Seek #Fno, 1 ReDim Qa(n + 1000) For I = 1 To n ReadRecord Fno: TurnRecord Nian, AX, 4 + 1, 1 n1 = Nian(1): m1 = Int(AX(1)): P = AX(2): E = AX(3): qb =AX(4)

Pe = P - k * E If Pe > 0 Then Eu = k * E: El = 0 XuManChanLiu1 Pe, R, Wm, Wmm, Wum, Wlm, W, Wu, Wl, b, X El R = 0 EuEl P, k, E, Wlm, Wu, Wl, Eu, El End IfRsRssRg1 Pe, R, S, Sm, Smm, Rs, Rss, Rg, Bx, XX, FR0, Kss,Kg, Immp EE = Eu + El: W = Wu + Wl HeWangZongRuLiu Rs, Rss, Rg, QRss0, QRg0, Qz, CRg, CRss, Qo,

U Tao = 1Cs = 1 - Cr * Qz ^ 0.4 Qa(I + Tao) = Cs * Qa(I) + (1 - Cs) * QzQ = Qa(I + Tao)Write #6 I, P, EE, R, Rs, Rss, Rg, W, Q, qb End IfNext IClo (6) = CStr(666)End Function

Public Sp(20) As VariantFunction ReadRecord(Fno) As Integer Dim TextLien As String, Stext As String, sss As String Dim I As Integer, j As Integer, LineLen As Integer j = 1

Line Input #Fno, TextLien LineLen = Len(TextLien) For I = 1 To LineLen Stext = Mid(TextLien, I, 1) Select Ca Asc(Stext) Ca 45, 46, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57 sss = sss + Stext If I = LineLen Then Sp(j) = sss End If Ca El If sss <> "" Then Sp(j) = sss: j = j + 1: sss = "" End If End Select Next I ReadRecord = jEnd Function

'蓄滿產流函數一Function XuManChanLiu1(Pe, R, Wm, Wmm, Wum, Wlm, W, Wu, Wl, b,X) As Integer Dim Dr As Single, A As Single A = Wmm * (1 - (1 - W / Wm) ^ X) If Pe + A < Wmm Then R = Pe - Wm + W + Wm * (1 - (Pe + A) / Wmm) ^ (1 + b) El R = Pe + W - Wm

End If If Wu + Pe - R < Wum Then Wu = Wu + Pe - R El Wl = Wl + Wu + Pe - R - Wum Wu = Wum If Wl > Wlm Then Dr = Wl - Wlm Wl = Wlm R = R + Dr End If End IfEnd Function

'兩層蒸發函數Function EuEl(P, k, E, Wlm, Wu, Wl, Eu, El) As Integer If Wu + P > k * E Then Eu = k * E Wu = Wu + P - Eu El = 0 El Eu = Wu + P El = (k * E - Eu) * Wl / Wlm Wu = 0 Wl = Wl - El End IfEnd Function

'三水源劃分函數Function RsRssRg1(Pe, R, S, Sm, Smm, Rs, Rss, Rg, Bx, XX, FR0,Kss, Kg, Immp) As Integer Dim j As Integer, nn As Integer Dim FR As Single, Au As Single, Smmf As Single, Smf AsSingle, RR As Single Dim Rs0 As Single, Rss0 As Single, Rg0 As Single, Kss0 As

Single, Kg0 As Single Rs = 0: Rss = 0: Rg = 0 If R > 0 Then FR = R / Pe S = S * FR0 / FR If FR < 1 Then Smmf = Smm * (1 - (1 - FR) ^ (1 / Bx)) El Smmf = Smm End If Smf = Smmf / (1 + Bx) nn = 5 RR = Pe / nn Kss0 = (1 - (1 - Kss - Kg) ^ (1 / nn)) / (1 + Kg / Kss) Kg0 = Kss0 * Kg / Kss For j = 1 To nn If S < Smf Then Au = Smmf * (1 - (1 - S / Smf) ^ XX) El Au = Smmf S = Smf

End If If (RR + Au > Smmf) Then Rs0 = (RR + S - Smf): S = Smf El Rs0 = (RR - Smf + S + Smf * (1 - (RR + Au) / Smmf)^ (1 + Bx)) S = S + RR - Rs0 End If Rs = Rs + Rs0 Rss0 = S * Kss0: Rss = Rss + Rss0 Rg0 = S * Kg0: Rg = Rg + Rg0 S = S - (Rss0 + Rg0) Next j FR0 = FR El Rs = 0 Rss = S * Kss Rg = S * Kg S = S - (Rss + Rg) End If If(Pe<=0) Then Pe=0Rs = Rs * FR0 * (1 - Immp) + Pe * Immp Rss = Rss * FR0 * (1 - Immp) ‘當Rs＝0時，自由水水庫面積仍前一時段的面積。 Rg = Rg * FR0 * (1 - Immp)End Function

'三水源河網總入流函數（地面徑流直接入流）

Function HeWangZongRuLiu(Rs, Rss, Rg, QRss0, QRg0, Qz, CRg,CRss, Qo, U) As Integer Dim Qs As Single, QRss As Single, QRg As Single Qs = Rs * U QRss = CRss * QRss0 + (1 - CRss) * Rss * U QRg = CRg * QRg0 + (1 - CRg) * Rg * U Qz = Qs + QRss + QRg + Qo QRss0 = QRss: QRg0 = QRgEnd Function