行政院國家科學委員會補助專題研究計畫
■ 成 果 報 告
□期中進度報告
集水區土地利用變遷對河川水文的影響
The Impact of Land Use Change on Stream Hydrology
計畫類別：□ 個別型計畫 ■ 整合型計畫
計畫編號：NSC 92－2621－Z －002 －005－
執行期間：
92 年 8 月 1 日至 93 年 7 月 31 日
計畫主持人：徐美玲
共同主持人：
計畫參與人員：
黃誌川、林庠序、杜佩儒、江瑩、徐慶安
成果報告類型(依經費核定清單規定繳交)：□精簡報告 ■完整報告
本成果報告包括以下應繳交之附件：
□赴國外出差或研習心得報告一份
□赴大陸地區出差或研習心得報告一份
□出席國際學術會議心得報告及發表之論文各一份
□國際合作研究計畫國外研究報告書一份
處理方式：除產學合作研究計畫、提升產業技術及人才培育研究計畫、列管
計畫及下列情形者外，得立即公開查詢
□涉及專利或其他智慧財產權，□一年■二年後可公開查詢
執行單位：國立台灣大學地理環境資源學系
中 華 民 國 93 年 10 月 25 日
中文摘要
近年來人口的增加和工、商產業的發展，使得都市邊緣可達性較高的丘陵山坡地逐漸 轉變為茶園、果園以及住宅等，這種土地利用的轉變對生態環境的衝擊，乃集水區經營者 在進行土地規劃和管理時必須面對的考驗。由於集水區內土地利用與自然環境因素之間的 交互作用，根據過去降雨和當時集水區土地利用狀況下所產生的逕流資料進行率定所得的 經驗型水文模式，無法有效推估集水區土地利用變動後的水文反應。因此本研究擬發展一 個適於本土環境的網格式分佈型水文模式（grid-based distributed hydrological model），以期 有效預測集水區土地利用改變時，所可能導致的河川水文特性的變化。目前已經完成模式中各相關模組的研發，包含利用DEMs 析取坡度、流向與集流面積；利用衛星影像評估土
地利用類別；結合蒙地卡羅法、最大概似法與效率係數的參數率定器；以及整合Green and
Ampt 的入滲理論、Manning 地表流速公式，Darcy 定律和擴散波理論的水文演算模組。
本年度承前兩年的研究結果，選擇面積為去年研究區7 倍的陳有蘭溪內茅埔上游集水 區作為研究區，從水利署所提供的近十多年時雨量資料中，選取五場降雨事件進行參數的 率定，並以二場事件進行驗證。由於缺乏實際觀測數據可以進行參數的校正，加上集水區 內的林地面積佔了九成多，因此將全區視為均質區域進行模擬。因為集水區面積增大，為 將模擬時間控制在合理的範圍內，而降低地形的空間解析度，加上降雨資料空間解析度的 不足，以及雨量紀錄經常出現缺損，致使無法有效掌握降雨的時空分佈特性，影響流量預 測的準確度。因此本模式雖對小集水區的水文模擬可以達到不錯的結果，但是當面積加大 空間變異也跟著增大時，雖能有效地模擬出主要的洪峰型態與大致的時間，但是流量變動 的誤差可能致使其實際利用價值降低。 關鍵字：地理資訊系統、數值地形模型、分散型水文模式、擴散波模式、粗糙係數、Green and Ampt 入滲理論
ABSTRACT
Rapid population growth and industrial development in recent years have resulted in urban sprawl onto the peripheral slopeland and cultivation of some forested land. These land use changes have great impacts on the watershed ecosystem, and have been the major concern in watershed management. Due to the interrelationship between the use of a land and its environmental factors, effective prediction of the impact of land use change on watershed runoff can only be achieved by incorporating the actual spatial distribution of land use within a watershed. This study aims to develop a grid-based distributed hydrological model which can readily incorporate land use information derived from satellite imageries and simulate watershed hydrographs for different land use spatial patterns. The model adopts continuity equation, diffuse wave model, and the Manning’s equation for surface runoff routing, while uses the Green and Ampt infiltration theory and Darcy’s Law for subsurface runoff treatment. Uncertainty analysis incorporating Monte Carlo simulation and a Maximum Likelihood Estimation and Efficiency Coefficient is used to retrieve optimal parameter values and to estimate the confidence interval for the simulation results.
The environmental data of the upper watershed of Chen-Yu-Lang River, and 5 rainfall events were used first for parameter calibration, and then 2 other storms were used to verify the results. Since the size of the study area is 6 times larger than that of the previous study area, it become necessary to lower the DEM resolution in order to facilitate the simulation at the expanse of topographic details. And as the study area increases, the variation of each hydrologic parameter increases accordingly, furthermore, inadequate rainfall data also result in poor simulation results. It is concluded that the proposed model is adequate in estimating hydrological responses in small watershed, however has shown its inadequacy for large watersheds.
Keywords: Geographic Information System, Digital Elevation Models, Distributed Hydrological Model, diffusive wave, Manning’s equation, Green and Ampt Infiltration.
目錄
目錄……… I
表目錄……… II
圖目錄………III
第一章 前言……… 1
第二章 研究方法……… 2
1. 模式架構 ……… 2
2. 研究步驟 ……… 3
（1） 研究區的選定 ……… 3
（2） 土地利用調查資料的蒐集………
3
（3） 降雨和水文資料的蒐集 ……… 3
（4） 相關地形參數的演算 ……… 4
（5） 入滲與地表下逕流 ……… 4
（6） 地表水流 ……… 4
（7） 模式的建構 ……… 6
（8） 模式的演算 ……… 6
（9） 模式參數的率定 ……… 9
（10） 模式的驗證 ……… 11
第三章 研究成果……… 12
1. 研究區概述 ……… 12
2. 水文資料 ……… 13
3. 土地利用與數值地形資料處理 ……… 14
4. 結果與討論 ……… 16
第四章 結論與建議……… 28
表目錄
表 1 研究區主要形貌參數……… 13
表
2 水文事件一覽表……….. 13
表 3 研究區各類土地利用面積與百分比……… 14
表
4 各作用參數率定值域與基本統計值………..17
表
5 各事件最佳參數組與綜合效率評估值………..17
表 6 率定與驗證事件的模擬結果評估表……….26
圖目錄
圖
1、模式架構與相關元件示意圖……… 3
圖
2 無限流向法網格分流示意圖 ……… 5
圖
3、不確定性分析流程示意圖……… 10
圖 4 研究區水系、地形、雨量站與徐昇氏分區圖……… 12
圖 5 研究區土地利用分佈圖 ……… 14
圖 6 研究區的坡度分佈圖……… 15
圖
7 研究區地形指數的空間分布圖……… 15
圖
8 1988 年 5 月 22 日降雨事件模擬結果與實際流量之比較……… 18
圖
9 1996 年 5 月 6 日降雨事件模擬結果與實際流量之比較………… 18
圖
10 1996 年 7 月 31 日降雨事件模擬結果與實際流量之比較………… 19
圖
11 1999 年 8 月 6 日降雨事件模擬結果與實際流量之比較 ………… 19
圖
12 2000 年 10 月 31 日降雨事件模擬結果與實際流量之比較 ……… 20
圖
13 2000 年 10 月 31 日降雨事件（未分組）模擬結果與實際流量之比較
圖
14 1998 年 8 月 3 日降雨事件模擬結果與實際流量之比較 ………… 21
圖
15 2000 年 8 月 22 日降雨事件模擬結果與實際流量之比較………… 21
圖
16 1998 年 8 月 3 日水文事件最佳參數組模擬地表水深時空分佈圖… 22
圖
17 1998 年 8 月 3 日水文事件最適參數組模擬地表水深時空分佈圖… 23
圖
18 1998 年 8 月 3 日水文事件最佳參數組模擬土體飽和度的時空分佈圖 24
圖
19 1998 年 8 月 3 日水文事件最適參數組模擬土體飽和度的時空分佈圖 25
圖
20 入滲容量在供水無虞條件與實際降雨變化比較圖……… 26
第一章 前言
許多水土災害的發生誘因皆與地表逕流或地表下逕流有關，如：土壤沖蝕與地表逕流、 淺層的崩塌與蠕移與地表下逕流皆有密切的關係。因此有效評估集水區內的地表與地下逕 流量及其分布將有助於其他相關研究的進行。以往通常是將集水區視為一個均質的空間單 元，然後利用相關理論與模式推估集水區的洪水流量，比較常用的包括：合理化公式、單 位歷線理論、瞬時地形單位歷線、SCS 的曲線數目法等；利用長期觀測資料的校定，這類 方法確能相當準確地推測集水區的洪水量，但是當集水區內部的地表特徵（土地利用的空 間分佈）發生變動時，卻往往導致推測值的誤差。40∼70 年代雖進行許多小型的試驗集水 區的研究，歸納出不同的土地覆蓋對於逕流量的影響。然而，試驗集水區的尺度較小(公頃 級)，往往無法直接評估中、大型集水區內土地利用變化所可能造成的影響。因此在 80 年 代初期，分散-物理型水文模式逐漸抬頭，將集水區拆解成許多小空間單元，根據地理資訊 系統提供的高程與土壤等空間資訊與水體的流動特性，計算集水區內部水文反應的空間分 布，如：SHE(Jonch-Clausen，1979；Beven et al.，1980)、SWAM(De Coursey，1982)、 IHDM(Morris，1979)等。這些模式大都沿襲著 Horton 的傳統，處理超滲地表逕流 (e.g. Huggins and Monke，1968；Smith and Woolhiser，1971；Engman and Rogowski，1974；Rovey et al.， 1977；Ross et al.，1979；Kutchment，1980)。然而，此類純物理型的水文模式對於 環境資訊（如土壤孔隙率、導水度等）的空間解析度要求較高，往往需要仰賴詳實的田野 調查，因此多僅適用於中小型的集水區。如今，隨著衛星科技(remote sensing)的發達與地 理資訊系統(geographic information system)的成熟，許多即時性的高解析度環境資料庫相繼 建立，使得進一步建立此等資料庫的資訊與水文參數之間的關係，以利經營管理單位有效 推估集水區水文特性的工作，更形重要（Beven，2000）。 為了避免大型集水區因局部大氣系統所導致的複雜性，如：雷雨胞移動方向與速度的 影響（徐義人，1995），及有效展現集水區內部空間變異特性，本研究乃針對中型規模的集 水區進行分散-物理型水文模式的建構。為了充分整合環境資料庫與水文模式，首先利用地 理資訊系統從現有的國土調查資料抽取土地利用的圖層，並轉換成相關水文參數，利用數 值地形（Digital Elevation Models，DEMs）析取相關的地形特徵，再以分散式水流演算 （Distributed Flow Routing）進行水流的演算。由於水文參數的空間變異性很大，如何有效 選擇代表性的參數一向是此類模式成敗的癥結，因此本研究演算時以蒙地卡羅模擬（Monte Carlo Simulation）選取相關水文參數，再進行模式不確定性的推算。一方面希望能將集水 區的土地利用和地形特徵有效地轉化成相關的水文參數，並盡量以物理模式進行水流的演 算，另一方面則希望能有效評估模式模擬結果的不確定性，以真實反應水文現象的高度空 間變異特性。第二章 研究方法
本研究擬以四年的時間進行分散-物理型水文模式的建構、測試與應用。第一年主要 係針對模式的各元件進行評估、測試及水文演算的模式原型建立。第二年著重在各元件的 整合及整體模式的實際測試以及相關誤差分析；第三年至第四年為模式的應用擬採用台灣 不同地區的不同尺度的集水區進行模式的應用。 1. 模式架構本研究根據數值地形模型（DEMs，Digital Elevation Models）計算影響集水區水流的 因素：集水區邊界、坡度、流向、地形指數，並利用地理資訊系統（GIS，Geographic Information System）整合土地利用資料，評估各種土地利用的可能入滲率與地表粗糙度，配合水流演 算（flow routing）建構推演分散型單場降雨事件的水文模式。利用 DEM 可以完整地描述 集水區內部高程(elevation)、坡度(slope)和坡向(flow direction)的空間變化，劃定集水區邊界 並推算各網格點的地形濕度指數(topographic wetness index)；根據 GIS 空間資料庫所提供的 土地利用資料，推估各網點的入滲特性（inf_a, inf_b）及地表粗糙度(na, nb)等四個與土地 利用相關的參數；藉由 Green and Ampt 的入滲理論估算入滲量，並以 Darcy-Weisbach equation 推估流速，再以聖凡南方程式確保網格間的水流連續性，並模擬集水區出口的水 文歷線，並利用蒙地卡羅模擬和效率係數組成的概似函數評估每一參數組模擬結果與觀測 值的概似程度，藉由貝氏定理獲得參數組的後驗機率分布，以趨近理想的參數數值，並利 用校準後的參數進行模式的驗證。研究架構如圖1 所示。 本研究第一年已完成根據數值地形模型演算相關地形參數、以地理資訊系統抽取土地 利用資料、利用相關參數進行地表水流演算等模式的初步建構。第二年加入地表下水流的 演算，並針對相關參數的率定方法進行研究，完成土地利用資料、颱洪事件的降雨量與逕 流量資料的收集與整理、模式效率係數的評估函式等的選取、參數率定方法的確立，加入 隨機變數的生成和模擬，並利用效率係數的計算，以實際降雨逕流資料進行模式的校正與 驗證工作。本年度則檢測前述完成模式運用在較大面積的集水區的可行性。
圖1、模式架構與相關元件示意圖 2. 研究步驟 （1）研究區的選定 為了進行模式有效性的驗證工作，必須選擇具有良好降雨和水文資料的集水區進行 研究；配合整合計畫的研究區，本年度選定南投縣內茅埔流量站以上的集水區為模式的 試驗集水區。 （2）土地利用調查資料的蒐集 為了獲得不同時期的土地利用資料，以探討研究區土地利用的變遷，擬以內政部地 政司的土地利用調查資料為基礎，蒐集研究區內的土地利用調查資料並加以彙整。該研 究區較完整的平地土地利用資料共有三期，分別為民國71 年、民國 77 年及民國 83 年的 全國土地利用調查資料；而山地地區僅有民國83 年的森林大調查一次資料。因此結合民 國 83 年的全國土地利用調查與森林大調查資料得到該集水區該年度土地利用的空間分 佈。 （3） 降雨和水文資料的蒐集 選擇研究區內及鄰近雨量站及流量站，自水利處官方網站收集並彙整該等測站的相 關降雨資料；利用徐昇式多邊形法進行集水區內部各網點雨量的空間內差。 降雨量時、空資料 數 值 地 形 土 地 利 用 模式模擬 二維漫地流，採用有限 差 分 法 之 二 次 交 互 疊 代 模 擬 每 一 時 間 之 漫 地流空間分布 模式校正 觀測事件 1. 水文歷線模擬 2. 土地利用變遷衝擊 坡度 流向 地形 指數 模擬 邊界 粗糙度評估 可滲性評估 空間內差 臨前評估 基流量、前期降雨
（4）相關地形參數的演算
Speight（1974）整理相關文獻發現，至少有 20 個以上與水文相關的地形參數可以利
用 DEM 資料估算出來，其中與分散型水文模式最相關者為坡度（slope）、坡向（flow
direction）與集流面積（contributing area）。本研究即利用台灣現有 40 公尺解析度的 DEM 計算各網格的坡度、坡向與集流面積。
（5）入滲（Infiltration）與地表下逕流（Subsurface flow）
在眾多的入滲理論中以Green and Ampt（1911）入滲理論較具物理意義，其以活塞流 (piston flow)的形式處理入滲水流向下移動的現象，也較易於操作，因此本研究採用此模 式，進行雨水的入滲演算： a inf P b inf t t b inf a inf f _ _ ； _ _ 0 − = + = （公式 1） 其中，f 為入滲率 [L/T]；inf_a、inf_b 兩參數乃受地表組成特性的影響，前者相當於
Horton 入滲理論中的穩定入滲率，後者則為控制入滲衰減率的參數；t0為ponding time[T]；
P 為降雨強度[L/T]。 雨水入滲所產生的地表下逕流則根據水力勢能來移動，以Darcy’s Law 來估計： dx dh K u sub sub = ⋅ （公式 2） 其中， K 為飽和水力傳導度[L/T]； hsub為水頭[L]；x 為水流距離[L]。 （6）地表水流 本研究以連續方程式和動量方程式進行網格間水流運動的推算，藉由 Green and Ampt 的入滲理論計算入滲量，利用曼寧公式估計地表水流速，並以聖凡南方程式來推算 各時間點網格間的水流。演算時採用的相關理論如下： 無限流向法 在網格式的資料結構下，計算水流方向及集流面積的方法大致上可分成單流向法和多 流向法兩類。單流向法最先是由 O’Callaghan and Mark（1984）所提出的，又稱為 D8 法 (deterministic eight-neighbors)，也是目前最簡單的計算方法。其乃以 3 ╳ 3 的移動視窗中， 鄰近8 格的最低點為中心網點水流的流向。後來 Quinn et. al.（1991）針對單流向失真的流 線缺失提出多流向法，其假設中心網格的水會流向周圍所有比它低的網格，而每個網格所 接收的水量則隨著其與中心網格的坡度而異(黃誌川等，2002)，各網格水流的分配可以下 式估計：
∑
= j j i i i L S L S f (公式 3)其中， S是八方位坡向的坡度，L是有效等高線長度（effective contour length），作為
加權之用；在東南西北的正交方向，L值取0.5，在對角鄰點則取 0.354。不過 L 值的大小，
不同地區的學者有不同的意見(Wolock et al.，1995；Freeman，1991；徐美玲，1995)。 目前利用上述兩種分流方法進行流量判定或集流面積推估的適用性尚未有明確的定論 （Lea，1992；Tarboton et al.，1997）。不過 Mattz and Garbrecht（1992）建議小尺度的邊坡 水文研究宜採用多流向法，但水流進入河道後則較適用單流向法。Tarboton（1997）整合
單流向與多流向的優點，提出無限流向法（圖2），經證實使用方便且結果較為正確（Zhou
and Liu，2002），因此本研究採用無限流向法進行水流的演算。
圖2 無限流向法網格分流示意圖 （改繪自 Tarboton, 1997） 擴散波模式（diffuse wave model）
任一時間點上，各點網格的水流必須符合連續方程式與動量方程式，水流速度的快慢 決定於各網格的勢能。Saint-Venant(1871)所提的明渠水流的動量方程式（公式 3）已成功 地用在漫地流的演算： 0 ) ( ) ( 1 1 0 2 = − − ∂ ∂ + ∂ ∂ + ∂ ∂ f S S g x y g A Q x A t Q A (公式 4) 其中，A：河槽斷面面積[L2]；Q：明渠流量[L3/T]；So：河床坡度；Sf：摩擦坡度；T： 時間[T]；x：水流方向[L]；y：垂直水流方向[L]；g：重力加速度[L/T2]
當水深對流速的影響可以忽略不計時，僅以河床坡度估算水流勢能的運算模式稱為運 動波模式（kinematic wave model），此種演算適用於坡度陡且回水效應不明顯的地形區。
若地形陡緩不一，則需考慮（公式 3）中的後兩項，即實際水面坡降及水流摩擦力，以獲
得較準確的水流估算，稱為擴散波模式(diffusion wave model)。由於研究區涵蓋坡地和河道 網格，因此採用此種演算法： x y S Sf o ∂ ∂ − = （公式 5） 有關水流摩擦力對水流的影響大多歸納自明渠實驗結果，其中 Manning’s equation、 Chezy equation 及 Darcy-Weisbach equaiton 最廣為採納，三者所採用的參數和原理大致相同 （Abrahams et al., 1991a, 1991b）。本研究則以曼寧公式進行水流流速的估算。
Manning’s Equation： N S h u 2 / 1 3 / 2 ⋅ = （公式 6） 式中，N為曼寧摩擦係數，主要受地面特性及水流深度的影響，可視為土地利用的函 數。 （7）模式的建構
以Green and Ampt 入滲理論分開入滲與超滲雨水，地表下逕流利用 Darcy’s Law，無 限流向法和連續方程式進行水流的演算，地表逕流則利用擴散波方程式、曼寧公式、無限 流向法和連續方程式進行各網格水流的演算： Continuity equation： f r x uh t h = − ∂ ∂ + ∂ ∂ ( ) （公式 7） 式中，t為時間[T]；h為格點的水流深度[L]；x 為水流方向的距離[L]；u 是該格點的 水流速度[L/T]；r為降雨強度[L/T]；f為格點的入滲率[L/T]。 （8）模式的演算---二階有限差分法
模式建構完成後，引入有限差分法（finite difference method）以DEMs的網格為架構， 解析上述各式，則地表逕流的連續方程式可寫成： ) ( 1 2 Qin Qout t h − δ = ∂ ∂ （公式8） Qin可寫為：
) ( Q Q 2 r f up in = +δ − （公式9） 地表下逕流則為： ) ( 1 2 Qin Qout t = − ∂ ∂ δ θ （公式10） Qin可寫為： f Qup in 2 Q = +δ （公式11） 式中的δ 是網格邊長[L]；Qin為該網格的總入流量[L3/T]；Qup為鄰近上坡網格的入流 量[L3/T]；Qout該網格總出流量[L3/T]。 本模式為提高顯性有限差分法的準確度，採用MacCormack（1969）所提出的兩步驟 修正預測法（two-step finite difference predictor corrector method）；以一次向前差分(Predictor)
與一次向後差分(Corrector)的平均求解。上述的各式改寫成： Predictor： 地表逕流收支平衡
(
1/2 1/2)
2 , , = + ⋅ + − + m i m i m i up m i in Q r f Q δ （公式12）(
m)
i out m i in m i i Q Q t h h* = + ∆2 , − , δ （公式13） 渠道流演算(
)
δ * * * ( ) , j j i i h Ele h Ele Sfi + − + = （公式14）( )
* 1/2 , 3 / 2 * * * 2 1 i f j i i i S R R N U ⋅ + ⋅ = （公式15） 漫地流演算 δ j i Ele Ele Sfi − = * , （公式16）( )
* 1/2 , 3 / 2 * * * 2 1 i f j i i i S h h N U ⋅ + ⋅ = （公式17） 地表下逕流收支平衡： 2 / 1 2 , _ , _ = + ⋅ + m i m i up sub m i in sub Q f Q δ （公式14）(
m)
i out sub m i in sub m i i Q Q t , _ , _ 2 * = + ∆ − δ θ θ （公式15）地表下逕流演算：
(
)
δ * , * , * ,i ( i subi) j subjsub Ele h Ele h
x h + − + = ∂ ∂ （公式16） x h K U subi i i sub ∂ ∂ ⋅ = * , * , （公式17） Corrector： 地表逕流收支： 2 * * * * * * ,i i i i i k out h h U A U Q = ⋅ = ⋅δ + （公式18） ) ( 2 1 * , , 2 / 1 , outm i outi m i out Q Q Q + = + （公式19） ) ( 1/2 1/2 2 2 / 1 , 2 / 1 ,+ =
∑
+ +δ + − + m i m i m j up m i in Q r f Q （公式20） ) ( 1/2 , 2 / 1 , 2 * * = + ∆ + − m+ i out m i in m i i Q Q t h h δ （公式21） 渠道流演算： ) ( 2 1 * ** 1 i i m i h h h + = + （公式22）( )
1 1/2 , 3 / 2 1 1 1 2 1 + + + + + = m i f m j m i i m i S R R N U （公式23） 2 1 1 1 1 1 1 , + + + + + + = ⋅ = m ⋅ im + km i m i m i m i out h h U A U Q δ （公式24） 漫地流演算： ) ( 2 1 * ** 1 i i m i h h h + = + （公式25）( )
1 1/2 , 3 / 2 1 1 1 2 1 + + + + + = m i f m j m i i m i S h h N U （公式26） 2 1 1 1 1 1 1 , + + + + + + = ⋅ = m ⋅ im + km i m i m i m i out h h U A U Q δ （公式27） 地表下逕流收支： 2 * , * , * , * , _ k sub i sub i sub i out sub h h U Q = ⋅δ + （公式28） ) ( 2 1 * , _ , _ 2 / 1 ,_ subm outi sub outi
m i out sub Q Q Q + = + （公式29） ) ( 1/2 , _ 2 / 1 , _ 2 * * = + ∆ + − m+ i out sub m i in sub m i i Q Q t δ θ θ （公式30）
地表下逕流演算： ) ( 2 1 * ** 1 i i m i θ θ θ + = + （公式31）
(
)
δ 1 , 1 , 1 ( ) + + + = ⋅ i + subm i − j + subm j i m i h Ele h Ele K U （公式32） 2 1 , 1 , 1 , 1 , _ + + + + = ⋅ + subm k m i sub m i sub m i out sub h h U Q δ （公式33） （公式12）至（公式33）中的m代表時間軸；∆t代表時間間隔[T]；為避免兩步驟交 互疊代過程中，因網格大小失當造成演算錯誤，另對時間間距及網格間距加以設限（Weinmann et al.，1979；Huang et al.，1985；Gharangik et al.，1991），其限制條件如下：
gh u t + ≤ δ ∆ 1 （公式34） 另外，就模式的初始和邊界條件的設定而言，在非河道的邊坡上，初始條件可視為無 水流狀態，即水深為零。台灣山地集水區由於坡度較陡，回水效應弱，因此集水區的流出 口假設上下網格水深一致（zero depth gradient condition），訂定其邊界條件（Morris，1979；
Govindaraju et al.，1990；Tayfur et al.，1993）。
（9）模式參數的率定
本模式首先由DEM擷取相關的地形參數，以並所蒐集的降雨─逕流模式中的水文參
數資料，做為本研究各相關參數的參考值，假設各參數在相對應的參考值範圍內呈分立 均等機率分布。然後利用由蒙地卡羅模擬（Monte Carlo Simulation）、最大概似法
（Maximum Likelihood Estimation）與改良的效率係數（Efficiency coefficient）組合而成
的參數率定方法－概似不確定性估計（Generalized Likelihood Uncertainty Estimation ，
GLUE；Beven，1992），進行模式中各類土地利用入滲（inf_a、inf_b）和粗糙係數（fri_a、
fri_b）等四個參數的選擇（圖3）。 根據各參數的範圍利用隨機變數生產器產生 5,000 組參數組合，以其分別進行模擬， 模擬所得流量的時間序列再與實測流量進行比對。由於本年度的研究區面積較大，且涵蓋 四個雨量站的降雨資料，初步僅利用效率係數（公式35）去評估各組參數對各事件模擬成 效時，發現對於水文歷線中最主要的洪峰流量和洪峰延時的準確度的評估權重太低，因此 調整效率係數的評估公式，增加洪峰流量和延時準確度的評估比重；另外以同組參數對其 他事件進行相同的模擬和評估，計算該參數組的概似函數（公式 36、37）（Lamb et al., 2000）。然後以結果最理想（也就是效率係數最高的）的100組繪出點狀圖與頻率分佈圖以 析取參數值域。另外，也藉由此100參數組分別繪出模擬結果，然後以上、下限百分之五
作為標準作為模擬結果之信賴區間推估。 圖3、不確定性分析流程示意圖 （公式35） （公式36） 上式中，θi為第 i 組參數組合，Y 所觀測的水文事件，Wi為一權重，依照模擬程度的優劣 主觀判定之，為簡單計算起見，概似函數多半介於0∼1之間。當水文觀測事件逐漸累積（由 Y1、Y2、Y3…Yi）之後，可利用貝氏定理加以修正，修正法如下：
∑
= = ⋅ + + ⋅ + ⋅ − ⋅ = M i N i N o N e N o e o e N i L L C W W W C Y L 1 ,...., 1 2 , 2 , 2 2 , 2 2 , 2 2 1 , 2 1 , 1 ,...., 1 ) | ( ... exp 1 ) | ( θ σ σ σ σ σ σ θ （公式37） 參數1 參數2 參數n 參數之初始機率 分立均等分配 水文模式模擬 利用M 個參數組 分別進行模擬 實測與模擬結果 的比較 概似程度評估 事後修正機率C為一權重調整係數，可避免數值的急遽上升。如此可計算出每一組參數組合針對 某場水文事件而言的概似程度；藉以分析這些參數組的可能分佈與模擬結果的最適區間。 洪峰流量的準確度乃以實測值（Pobs）與模擬值（Psimu）的比值為代表。 % 100 × = obs simu p P P E 當 Psimu ≤Pobs % 100 × = simu obs p P P E 當 Pobs ≤Psimu 洪峰延時的準確度，則以實測洪峰到達時間與模擬時間相差3小時以上者為0，2小 時者為34％ ，1小時者為67％，而完全沒有差異者為100％。 將洪峰的流量的準確度、洪峰延時的準確度和整體水文歷線的準確度（EFFIC），分 別以0.5、0.3和0.2為權重值，相加得到綜合效率的評估數值。 （10）模式的驗證 利用選定的參數組進行其他水文事件的模擬，比對模擬的水文歷線與實際觀測情 況，針對洪峰流量的大小、遲滯時間和歷線的效率係數進行誤差的檢核，以評估模式的 有效性。
第三章 研究成果
1. 研究區概述 位於南投縣信義鄉境內的陳有蘭溪發源於玉山北坡，東以郡大山與郡大溪流域相接， 西以阿里山與清水溪及曾文溪流域為鄰，南以玉山、鹿林山與高屏溪流域分隔，匯集玉山 以北、阿里山山脈東坡與郡大山脈西坡的水，由南向北流。本研究區乃涵蓋水利署位於內 茅埔溪匯流口的內茅埔流量站以上的陳有蘭中上游集水區，總面積達368.43平方公里。區 內最高海拔達2900公尺，最低處約670公尺，標高1000至3000公尺的土地約佔研究區面 積之半，區內邊坡坡度以30至 40度者最多，約佔37％，而平均坡度約為34度。研究區 內主流河段長約27820公尺，平均坡降約為0.76。 研究區的地質大致以陳有蘭溪及其支流沙里仙溪為分界，東側為古第三紀變質岩區， 以硬頁岩、板岩、千枚岩、變質砂岩和石英岩為主；西側則為新第三紀沈積岩，以砂岩、 頁岩和砂頁岩互層為主。研究區包括沙里仙溪斷層、東埔斷層等主要的構造線，以及和社 背斜、同富山向斜、神木斷層等，構造線附近的岩體受到構造應力作用而特別破碎，風化 侵蝕作用也較為旺盛。 圖 4 研究區水系、地形、雨量站與徐昇氏分區圖 表 1 研究區主要形貌參數變數 面積 （km2） 周長 (km) 密集度 平均坡度 (%) 平均高程 (meter) 最遠流長 (km) 研究區 367.43 89.54 1.65 60.67 1715.92 35.18 2. 水文資料 研究區內每年5月到6月的梅雨季節的降雨量約佔全年的30％，而7月到9月的颱風 季節則降下全年雨量之半；因此5月到9月的雨季降雨量約佔全年降雨量的70％左右（陳 振杰，1999）。研究區出口處的隸屬水利署的內茅埔流量測站提供本研究所需的時流量資 料，另外，區內內茅埔、望鄉、東埔以及鄰近的新高口雨量測站提供本研究所需的時雨量 資料。 本研究將從水利署取得近十年來的時雨量與時流量資料，以5小時連續不降雨為切割 標準，取得7個總雨量較大且資料完整的水文事件，各事件的基本統計特性如表2所示： 表2 水文事件一覽表 日期 歷時 (hr) 總降雨量 (mm) 洪峰稽延 (hr) 洪峰流量 (cms) 總逕流量 (cms) 逕流係數 (%) *1988/05/22 59 161.38 3 680 9819.9 59.6 *1996/05/06 45 91.36 2 225 3169 34.0 *1996/07/31 43 813.77 0 3370 60673.0 73.1 *1999/08/06 58 87.73 4 298 1837.6 20.5 *2000/10/31 33 147.68 2 248 3895 25.8 **1998/08/03 60 144.01 2 642 7532 51.2 **2000/08/22 60 186.8 16 479 8970 47.1 逕流率僅供參考；*為模式率定事件；**為模式驗證事件。 由於本研究區規模約達 370平方公里，約為上期研究區的 3.5 倍，因此雨量的處理較 為複雜，為了比對洪峰與雨峰的時間差異，乃利用徐昇氏法將研究區劃分為四個次區（圖 4），以其面積作為權重加總內茅埔、望鄉、東埔以及新高口雨量測站的時雨量資料以獲得 一綜合雨量主體圖（圖 8 至圖 15），然後比較洪峰與綜合雨峰的時間差。結果 1995 年 7 月1日的降雨事件所得的洪峰發生時間正好與綜合雨峰一致；，而2000年8月22日因為 的逕流有兩個高峰，比較高的洪峰與主要雨峰之時間差達16小時，但是第一個洪峰則與綜 合雨峰差距不到一小時。另外對於某些水文事件而言，水利署所提供的流量資料持續時間 不夠，可能因而導致總逕流量的低估。又因為利用徐昇氏法切割所得各雨量測站所代表的 次集水區，可能並不能全然代表實際降雨的空間分佈，因此導致逕流係數的大幅變動。
3. 土地利用與數值地形資料處理
由於本研究乃選取發生在1988∼2000年間的水文事件，因此相對應的土地利用資料，
平地部分參考1994年內政部地政司所進行的五千分之一國土普查資料，山區部分則以1994
年農委會所進行的森林大調查資料為主。利用 Arc_View 將此兩份資料進行屬性表合併，
並依照Kaltsas and Bosch (2000)的分類標準區分為五大類（表3），以利模式中各參數的套
疊，並將所得向量圖層轉換為與數值地形檔相同大小的網格式資料。各類土地利用的空間 分佈與所佔的百分比如圖5及表3所示。因為林地約佔整個集水區的93％，且沒有兩次土 地利用資料可供率定不同土地利用之間的參數，因此本研究乃將全區視為均一的土地利 用，以利參數不確定性的分析。 表 3 研究區各類土地利用面積與百分比 土地利用 水體 open water 森林 forest 農作物 herbaceous 人為用地 distributed 混合地 mixed 面積(km2) 7.41 340.09 14.51 1.35 3.08 比例(%) 2.01 92.56 3.95 0.37 0.84 圖 5 研究區土地利用分佈圖 另外，利用數值地形資料析取坡度、流向、地形指數、集水區邊界、河道等五項物理 參數的基礎；其中，坡度採用三次加權法，所得坡度分布如圖6所示：
圖 6 研究區的坡度分佈圖 本研究採用無限流向法計算各網點的集流面積，再以Arc_View的map_calculator計算 地形指數，圖7顯示研究區的地形指數空間分佈。採用TOPMODEL的方法，以集水區內 所有格點的地形指數的平均值作為模式初始模擬時的土體飽和度判定準據，將數值在平均 值以上的格點視為飽和；而在平均值以下的格點則以其與平均值的比率為該點的飽和度。 圖7 研究區地形指數的空間分布圖
N
N
4. 結果與討論 本研究利用不確定分析方法來進行參數值域的析取，為減少模擬時所需的計算資源， 先針對參數值域的合理性進行評估，利用隨機變數生產器產生參數值域較為寬廣的500組 均等分佈的參數組合，分別為：飽和水力傳導度（K，m/min）介於 0.1∼0.000001、孔隙 率（n，%）介於 0∼1、飽和度的標準（S）介於 0∼30、粗糙係數（N）介於 0∼1、穩定 入滲率（In_a，mm/hr）介於0∼70、入滲變化率（In_b，mm）介於0∼70，利用研究集水 區 15 場水文事件進行模擬，再將所得結果的的分析，調整各參數值域為：K 介於 0.1∼ 0.001、n介於0∼0.8、S介於0∼25、N介於0∼0.2、In_a介於0∼70、In_b介於0∼70， 以做為隨機變數生產器的值域範圍。 選擇五場水文事件做為模式率定的基礎（表2），分別包含不同降雨規模、洪峰流量與 不同基流狀況。本研究區的DEMs總網格數為339*378 = 128,142個格點，模擬一場30小 時的水文事件大概約需 3.5 分鐘（電腦配備為：CPU = Pentium_4 1.4 GHz、Ram = 512
MB），而模擬一場45hr 的水文事件約需5 分鐘，模擬一場60 小時的水文事件大概約需7 分鐘。換句話說，完成五場事件各場5,000次模擬所需的時間約需103天。 利用初步選定的參數區間及分立均等分配，以亂數抽取生產 5,000 組新的參數組進行 模擬，然後計算各組參數模擬結果的效率係數。從效率係數與各參數值的關係圖中，可以 看出各水文事件在 5,000 次的模擬中都有一些表現不錯的參數組，並顯現出各參數的特性 與較為合理的值域。為合理歸納選取有效的參數組，從每場事件 5,000 組的模擬中，先選 取綜合效率評估值在70％以上者，得60組參數組，然後再從各組中評估值不到70％者挑 10組表現最佳者，篩選掉先前已經選取者後，再從各組挑5組，直到總參數組達100組。 此100組中各參數的值域範圍如表4所示。 本模式首先以這最佳的100組參數組分別針對單、雙雨峰的校正事件進行模擬，從中 各選平均表現最佳之一組，做為該組的最適參數組(單峰：K =0.028 m/hr、n = 0.61、S = 11.21、N = 0.0281、In_a = 58.36 mm/hr、In_b = 25.74 mm；雙峰：K =0.058 m/hr、n = 0.54、 S = 25.8、N = 0.0068、In_a = 11.36 mm/hr、In_b = 15.21 mm)。圖8到圖12呈現針對5場 校正水文事件利用此參數組及各場最佳參數組（表 5）所做的模擬水文歷線與實際觀測流 量的比較。其中1996年7月31日賀伯颱風的異常降雨量超過本模式的計算容量，因此出 現洪峰被截斷的情況（圖10）；而1996年5月6日的降雨事件有兩個明顯的洪峰，雖然第 一個洪峰流量較高，但是第二個洪峰維持較久（圖 9）。本研究非常明顯地掌握到兩次洪峰， 若以所得洪峰與第一個洪峰去評估降雨延時差異較大，但若與第二個洪峰則顯然非常準 確，使最佳效率評估值提高為 0.78（表 5）。至於其他的降雨事件最佳效率評估值都在 75 ％以上，甚至有達 98％者，應該算是不錯的結果。圖 13 顯示根據雨型分組可以提高最適 組參數的模擬效率評估值，不過從最適值均較最佳值為低的情況來看，取得均吻合各場降 雨事件的參數組並不容易。
表4 各作用參數率定值域與基本統計值 基本統計 參數 最適值 平均值 中位數 標準差 上限 下限 K (m/hr) 0.028 0.040 0.028 0.030 0.094 0.001 n (%) 0.61 0.46 0.51 0.27 0.89 0.06 S 11.21 13.82 12.58 7.61 26.89 1.14 N 0.0281 0.155 0.037 0.203 0.527 0.001 In_a (mm/hr) 58.36 28.33 20.82 20.47 60.79 4.82 In_b (mm) 25.74 35.34 32.86 20.87 66.61 3.63 表5 各事件最佳參數組與綜合效率評估值 綜合效率評估值 降雨事件 K n S N In_a In_b 最佳 最適 1988/05/22* 0.0428 0.63 10.90 0.0232 47.49 18.16 0.75 0.68 1996/05/06* 0.0278 0.61 11.21 0.0281 58.36 25.74 0.63（0.78） 0.63（0.78） 1996/07/31 0.0746 0.32 20.64 0.1056 36.67 62.07 0.47 0.45 1999/08/06 0.0602 0.94 11.27 0.0390 46.72 66.61 0.98 0.78 2000/10/31* 0.0693 0.57 20.16 0.0014 15.93 21.23 0.88 0.69 1998/08/03 0.0635 0.58 12.28 0.0124 54.69 33.21 0.90 0.72 2000/08/22* 0.0536 0.90 12.25 0.0051 46.82 42.18 0.50（0.65） 0.22（0.37） 前5場為校正事件，後兩場為驗證事件；＊表示該場降雨事件有兩個明顯的雨峰 利用前述選定的參數組，以剩下的2場水文事件進行前述選定參數組的驗證，並以最 佳100組對驗證降雨事件進行驗證，將所得的流量曲線中，共同時間點的模擬流量值，去 掉最大和最小各5個數值，以訂定該時間點模擬區間（圖14∼圖15）。並分別利用最適參 數組與最佳參數組進行模擬，得到研究區地表水深與土壤飽和度隨著時間變化的空間分佈 圖（圖16到圖19）。另將此等參數組對各校正事件與驗證事件的模擬結果與實際流量歷線 進行比對與評估（表6 ）。
圖8 1988年5月22日降雨事件模擬結果與實際流量之比較 圖9 1996年5月6日降雨事件模擬結果與實際流量之比較 0 5 10 15 20 25 30 35 40 45 50 0 5 10 15 20 25 30 35 40 45 50 55 60 time(hr) ra in fa ll( mm) 0 100 200 300 400 500 600 700 800 Q( c m s) 0 5 10 15 20 25 0 5 10 15 20 25 30 35 40 45 time(hr) ra in fa ll( mm) 0 50 100 150 200 250 Q( c m s)
圖10 1996年7月31日降雨事件模擬結果與實際流量之比較 圖11 1999年8月6日降雨事件模擬結果與實際流量之比較 0 10 20 30 40 50 60 70 80 0 5 10 15 20 25 30 35 40 45 50 55 60 time(hr) ra in fa ll( mm) 0 500 1000 1500 2000 2500 3000 3500 4000 Q( c m s) 0 5 10 15 20 25 0 3 6 9 12 15 18 21 24 27 30 time(hr) ra in fa ll( mm) 0 50 100 150 200 250 300 350 400 450 Q( c m s)
圖12 2000年10月31日降雨事件模擬結果與實際流量之比較 圖13 2000年10月31日降雨事件（未分組）模擬結果與實際流量之比較 0 5 10 15 20 25 0 5 10 15 20 25 30 35 40 45 50 55 60 time(hr) ra in fa ll (m m ) 0 200 400 600 800 1000 1200 Q( c m s) 0 5 10 15 20 25 0 5 10 15 20 25 30 35 40 45 50 55 60 time(hr) ra in fa ll( m m ) 0 200 400 600 800 1000 1200 Q( c m s)
圖14 1998年8月3日降雨事件模擬結果與實際流量之比較 圖15 2000年8月22日降雨事件模擬結果與實際流量之比較 0 10 20 30 40 50 60 0 5 10 15 20 25 30 35 40 45 50 55 60 time(hr) Ra in fa ll( m m ) 0 200 400 600 800 1000 1200 Q( cm s) 0 5 10 15 20 25 0 5 10 15 20 25 30 35 40 45 50 55 60 time(hr) R a in fa ll( m m ) 0 200 400 600 800 1000 1200 Q( c m s)
圖16 1998年8月3日水文事件最佳參數組模擬地表水深時空分佈圖 Time = 1 hr Time = 5 hr Time = 10 hr Time = 15 hr Time = 25 hr Time = 40 hr 0 4.5 9 18 27 36 Kilometers
N
圖17 1998年8月3日水文事件最適參數組模擬地表水深時空分佈圖 Time = 1 hr Time = 5 hr Time = 10 hr Time = 15 hr Time = 25 hr Time = 40 hr 0 4.5 9 18 27 36 Kilometers Time = 5 hr
N
圖18 1998年8月3日水文事件最佳參數組模擬土體飽和度的時空分佈圖 Time = 1 hr Time = 5 hr Time = 10 hr Time = 15 hr Time = 25 hr Time = 40 hr 0 4.5 9 18 27 36 Kilometers
N
圖19 1998年8月3日水文事件最適參數組模擬土體飽和度的時空分佈圖 Time = 1 hr Time = 5 hr Time = 10 hr Time = 15 hr Time = 25 hr Time = 40 hr 0 4.5 9 18 27 36 Kilometers
N
表 6 率定與驗證事件的模擬結果評估表 評估結果 水文事件 效率評估 值（％） 效率係數 (%) 洪峰時 間誤差 (hr） 實測洪 峰流量 (cms) 模擬洪 峰流量 (cms) 洪峰流 量誤差 (%) 1988/05/22* 64(68) 59 (-57) 1 680 546.4 -19.6 1996/05/06 78 57 8（0） 225 214.5 -4.7 1996/07/31 45 -14 0 3370 1072.3 -68.2 1999/08/06 78 33 0 298 428.2 43.7 2000/10/31* 69(42) -1681 (-7) 1 248 370.8 49.5 1998/08/03 72 22 0 642 1075.0 67.4 2000/08/22* 37(22) -340 (-300) 13（0） 479 739.9 54.5 ● 效率評估值與效率係數欄位中的數字表示利用所以校正事件所得最適參數組與分降雨組別後的最適組 （比括號表示）模擬結果。 ● 洪峰時間差中（）表示第二個洪峰時間的預測誤差 導致模擬誤差的原因不只一端，就模式而言，本模式利用曼寧粗糙係數和 Green and Ampt入滲參數去調節水流入滲和流動的速度。地表水流的流速乃由水力半徑（或水深）、
坡降與粗糙係數去推估，以往的實驗證實(Abrahams et al.，1991a；1991b)，粗糙係數與渠 道的形狀或當時的水深（或水力半徑）有關，但是本模式為求演算的可行性，每次演算時 將整個集水區的粗糙係數視為定值；換句話說，原本水深較高的格點，水流的實際流速應 該很大，讓降雨間歇時水很容易消退，然而因為模式的限制未能讓粗糙係數隨著水深調整， 因此可能因而導致水流速度過慢。
另外，Green and Ampt（1913）的入滲理論乃根據供水無虞的條件所做的運算與推導，
而自然降雨卻非常複雜，當降雨間歇時，入滲率應該不會依照原本的入滲曲線繼續下降， 很可能因土壤排水而再度增大，因此導致計算的誤差。而此種入滲率的變化可由 Hewlett (1967)的野外觀察中獲得證實（圖20），但由於此種變動尚未有具體的修正模式，因此本模 式目前也無法有效模擬降雨間歇時所產生的流量變動。 圖20 入滲容量在供水無虞條件與實際降雨變化比較圖 (Hewlett et al.，1967)
然而最大的問題癥結有二，一是集水區面積過大，一是雨量資料不足。當面積過大時， 坡地水流至河道時水量過大，而模式中的河道寬度僅限一格（80 公尺），當實際河道因洪 水上漲，致使河道中的水深快速上漲，甚至超過兩岸邊坡格中的水位時，不但該點坡地水 流無法匯入，而且由於模式在水流進入河道後採用單流向演算，無法往兩岸邊坡格點散流， 因此一旦河道流量積高太多，演算結果其實已經偏離實際情況，導致模式的中斷。又因為 本研究區面積過大，為進行模擬工作而將地形解析度由40公尺降低至80公尺，可能使得 地形失真。而區內各參數的變異隨著面積增大而更形提高，目前由於沒有小區域的觀測資 料，無法進行實際的參數校正，因此僅能以單一的數據進行模擬，導致的模擬誤差自然比 小面積的模擬誤差大得多。 本研究共使用四個雨量測站的資料，但是由於某些測站僅從1982年以後開始有時雨量 資料，加上許多測站的資料不全，使得可用的降雨事件降低。而運用徐昇氏武斷區劃各測 站雨量的代表範圍，不但面積過大，而且與地形分區並不一致，因此不能有效模擬研究區 降雨的時空分佈狀態。另外，全區各參數均以單一數值進行模擬，結果分別以各水文事件 的洪峰流量、延時和流量歷線所率定的最佳模擬參數組並不一致，因此最適模擬參數組的 選取受到各參數主觀評估權重的影響。
第四章
結論與建議
從理論剖析與實證研究的結果顯示，本研究所建構的水文模式可以有效模擬小集水區 的水文歷線、洪峰時間、洪峰流量，與地表下逕流，提供地表水深與土體飽和度隨時間變 化的空間分佈。但是由於各參數與土地利用之關係尚未釐清，仍僅能以全區均質的方式進 行模擬，而且當集水區面積增大時，為將模擬時間控制在合理的範圍內，必須降低地形的 空間解析度，加上降雨資料空間解析度的不足，使得模擬準確度大幅降低。就小集水區而 言，或許可以達到不錯的預測效果，但是當面積加大時空間變異也跟著增大，再以全區共 用參數去模擬流量可能就不盡理想，雖能有效地模擬出主要的洪峰型態與大致的時間，但 是流量變動的誤差可能致使其實際利用價值降低。未來有關的研究應該仍以小集水區較為 可行，聚焦在各類土地利用與水文參數之間關係的探討，確立效率評估參數的明確準則， 以及率定和建立模式使用於不同地理區時的合理參數資料庫。 另一個模擬結果的誤差乃導源於雨量資料的不足。台灣的地形崎嶇，降雨的時空分佈 變異甚大，但是目前山區的雨量站寥寥可數，根本無法有效顯示各場降雨事件的時空分佈 特性。而許多有關的災害，如土石流、洪水的預警又多仰賴有效降雨資訊的掌握，因此山 區雨量站的廣泛裝設應該是刻不容緩的當務之急。而如何有效將少數幾點的雨量時序資料 轉換成可以操作又合理的時空變動資訊，也是未來值得進一步探索的研究議題。第五章 引用文獻
徐美玲 (1995) “坡地土壤孔隙水壓動態空間分佈預測模式”，國立台灣大學地理學系地理學
報 18：1-22。
黃誌川、徐美玲、楊奕岑 2002 多流向與單流向集流面積估算結果之空間分佈特性，中華
水土保持學報，33 (1)：57-72。
Abrahams, A. D. and Parsons A. J., 1991a, “Relation Between Infiltration and Stone Cover on A Semiarid Hillslope”, Southern Arizona, Journal of Hydrology, 122: 49-59.
Abrahams, A. D. and Parsons, A. J., 1991b, “Resistance to Overland Flow on Desert Pavement and Its Implications for Sediment Transport Modeling”, Water Resources Research, 27(8): 1827-1836.
Beven, K. 2000. On the future of distributed modeling in hydrology, Hydrological Processes, 14: 3183-3184.
Beven, K. J., Warrwn, R. and Zaoui, J. 1980. SHE: towards a methodology for physically-based distributed forecasting in hydrology, In Hydrological Forecasting. IAHS Publication 129, 133-137.
DeCoursey, D. G. 1982. ‘ARS’ Small watershed model, Am. Soc. Agric. Eng. Paper, 2082-2094. Engman, E. T. and Rogowski, A. S. 1974. A partial area model for storm flow synthesis, Water
Resources Research, 10(3): 464-472.
Freeman, T. G., 1991, “Calculating catchment areas with divergent flow based on a regular grid”,
Computers and Geosciences, 17(3): 413-422.
Gharangik, A. M. and Chaudhry, M. H., 1991, “Numerical simulation of hydraulic jump”,
Journal of Hydrologic Engineering, 117(9): 1195-1232.
Govindaraju, R. S., Kavvas, M. L. and Jones, S. E., 1990, “Approximate analytical solutions for overland flow”, Am. Geophys. Union, Water Resources Research, 26(12): 2903-2912. Green, W. H. and Ampt, G. A., 1911, “Studies on soil physics, part I, the flow of air and water
through soils”, Journal Agriculture Science, 4(1): 1-24.
Hewlett, J. D. and Hibbert, A. R. 1967 Factors affecting the response of small watersheds to precipitation in humid areas. In: Inter. Symp. on Forest Hydrology, Sopper and Lull (Eds), Pergamon Press, Oxford, pp. 275-290.
Huang, J. and Song, C. C. S., 1985, “Stability of dynamic flood routing schemes”, Journal of
Hydrologic Engineering, ASCE, 111(12): 1497-1506.
Huggins, L. F. and Monke, E. J. 1968. A mathematic model for simulating the hydrologic response of a watershed, Water Resources Research, 4(3): 529-539.John—Clausen, T. 1979. Systeme Hydrologique Europeen: a short description. SHE Report I, Danish Hydraulics Institute, Horsholm, Denmark.
Jain, S. K., Singh, R. D. and Seth, S. M. 2000. Design flood estimation using GIS supported GIUH approach, Water Resources Management, 14: 369-376.
Kutchment, L. S. 1980. A two-dimensional rainfall-runoff model: identification of parameters and possible use for hydrological forecasts. In Hydrological Forecasting, IAHS Publication 129: 215-219.
Lea, N. J., 1992, “An aspect-driven kinematic routing algorithm, in overland flow”, In A. J. Parsons and A. D. Abrahams(eds.), Hydraulic and Erosion Mechanics, Chapman and Hall, New York, 1992.
MacCormack, R. W., 1969, “The effect of viscosity in hypervelocity impact cratering”, Am. Inst.
Aeronaut and Astronaut., N. Y., 69-354.
Martz, L. W. and Garbrecht, J., 1992, “Numerical Definition of drainage network and subcatchment areas from digital elevation models”, Computer and Geosciences, 18(6): 747-761.
Morris, E. M. 1979. The effect of the small-slope approximation and lower boundary conditions on solutions of the Saint-Venant equations, Journal of Hydrology, 40: 31-47.
O’Callaghan, J. F. and Mark, D. M. 1984. The extraction of drainage networks from digital elevation data, Graphic and Image Processing, 28: 323-344.
Quinn, P. F., Beven, K. J., Chevallier, P., and Planchon, O., 1991, “The prediction of hillslope flow paths for distributed hydrological modeling using digital terrain models”,
Hydrological Processes, 5: 59-79.
Ross, B. B., Contractor, D. N. and Shanholz, V. O. 1979. A finite element model of overland and channel flow for assessing fthe hydrological impact of land-use change, Journal of
Hydrology, 41: 11-30.
Rovey, E. W. and Woolhiser, D. A. 1977. A distributed kinematic model of upland watersheds,
Smith, R. E. and Woolhiser, D. A. 1971. Overland flow on a infiltrating surface, Water Resources
Research, 7(4): 899-913.
Speight, J. G.,. 1974, “A parametric approach to landform regions”, Special Publication Institute
of British Geographers, 7, 213-230.
Tarboton, D. G., 1997, “A New Method for the Determination of Flow Directions and Contributing Areas in Grid Digital Elevation Models”, Water Resources Research, 33(2): 309-319.
Tayfur, G., Kavvas, M. L., Govindaraju, R. S. and Storm, D. E., 1993, “Applicability of St. Venant equations for two-dimensional overland flows over rough infiltrating surfaces”,
Journal of Hydrologic Engineering, ASCE, 119(1): 51-63.
Weinmann, P. E. and Laurenson, E. M., 1979, “Approximate flood routing methods: A review”,
Journal of Hydrologic Engineering, ASCE, 105 (12): 1521-1537.
Wolock, D. M. and Mc Cabe Jr., G. J., 1995, “Comparison of single and multiple flow direction algorithm for computing topographic parameter in TOPMODEL”, Water Resources
Research, 31(5): 1315-1324.
Zhou, Q. and Liu, X., 2002, Error assessment of grid-based flow routing algorithms used in hydrological models, International Journal of Geographic Information Science, 16(8): 819-842.
