國立臺灣大學土木工程學研究所 碩士論文
Department of Civil Engineering College of Engineering
National Taiwan University Master Thesis
基隆河上游集水區含員山子分洪道之出水口水位 預報模式
Forecasting Outlet Water Level with Yuansantze Flood Diversion Tunnel for Keelung River Upstream Watershed
邱啟平 Chi-Ping Chiou
指導教授:蔡丁貴 教授 Advisor: Ting-Kuei Tsay, Prof.
中華民國 98 年 6 月
June, 2009
中 文 摘 要
由於基隆河上游興建員山子分洪道設施,將颱風或豪大雨造成的洪水進行部 分疏導,使得原有的水位預報模式(吳等,2004)計算大華橋水位站有高估水位的 現象。本研究目的,在於考量員山子分洪量等因素並修正模式,使大華橋水位能 夠更精準的預報,其預報結果可提供下游洪水預警,爭取沿岸居民疏散時間。模 式主要考慮分洪量、河川水位之變動與其臨前狀態之序率關係,同時加入可能影 響水位變動之因素:集水區之降雨量。本文假設河川水位在某一個時刻,為其所 有影響因子的線性組合,使用最小平方法,從歷史颱風或豪大雨事件之水文記錄 中率定出一條具有迴歸關係特性的預報函數式,其函數式再經由模擬退火法修 正。當豪大雨事件發佈時,蒐集其集水區內各站之降雨量、員山子水位等等資訊,
利用此一時序性的序率遞迴關係函數式即可預報該水位站一段時間之水位變化。
本模式應用柯羅莎(Krosa)、薔蜜(Jangmi)等颱風預報大華橋水位之變化,皆獲得 良好之驗證。
關鍵詞:序率、最小平方法、模擬退火法、水位預報
Abstract
Yuansantze Flood Diversion Tunnel is located at the upstream of Keelung River.
The Tunnel diverts part of flood from upstream basin discharge during typhoon periods.
It results in over-estimation of forecasted water levels at the Dahua Bridge station (Wu.
et al., 2004). This paper takes into accounts of flood diversion to extend capability of previous forecast model. The accurately predicted results at the Dahua Bridge Station can provide with better forecasting of downstream flood water levels. This improvement will allow inhabitants residing along river to be evacuated timely.
Stochastic relation includes flood diversion, present and the antecedent records at the river gage, and rainfalls are primarily concerned in the model. This paper assumes that water stage at the outlet of watershed at any time is linear combination of all the affecting factors. Based on historical data during typhoons events, a recursive relationship is developed by employing the least squares method and simulated annealing algorithms. With the recursive relationship formula, water levels at the Dahua Bridge Station can be predicted more accurately. Present model is applied to predict water levels at the Dahua Bridge Station for the Typhoon events of Krosa and Jangmi. Good agreement between forecasted and measured results is observed.
Keywords: Water level forecast, Stochastic, Least squares method, Simulated annealing algorithms
目 錄
中文摘要……… I 英文摘要……… II
第一章 緒論 ... 1
1.1 前言 ... 1
1.2 文獻回顧 ... 2
1.2.1 降雨逕流或水位推估 ... 2
1.2.2 模擬退火演算法 ... 2
1.3 研究目的 ... 3
第二章 模擬退火演算法之理論 ... 5
2.1 模擬退火法之原理及物理現象 ... 6
2.2 模擬退火法之搜尋方式及演算步驟 ... 8
2.3 模擬退火演算法影響因子之探討 ... 12
第三章 模式建立 ... 15
3.1 序率分析法 ... 15
3.2 最小二乘法 ... 17
3.3 模擬退火演算法於集水區參數調整之運用 ... 18
3.4 建置步驟 ... 21
第四章 模式驗證 ... 25
4.1 研究區域 ... 25
4.2 集水區出水口大華橋之水位預報模式驗證 ... 28
4.2.1 模式率定 ... 29
4.2.2 模式驗證 ... 32
4.3 淡水河全流域不恆定流河川水位預報模式驗證 ... 36
4.4 員山子分洪流量 ... 43
4.4.1 員山子分洪工程概述 ... 43
4.4.2 馬斯金更法(Muskingum method) ... 48
第五章 結論與建議 ... 54
5.1 結論 ... 54
5.2 建議 ... 56
參考文獻 ... 57
附錄一 全流域不恆定流模式介紹 ... 60
A.1 模式理論 ... 60
A.1.1 控制方程式 ... 60
A.1.2 第二種多方式特徵法 ... 62
A.2 邊界條件 ... 63
A.3 全流域不恆定流模式架構 ... 63
A.4 模組限制及注意事項 ... 68
A.5 模組應用 ... 72
表目錄
表2.1 物理系統與搜尋最佳解類比對照表 ... 8
表3.1 相關參數設定 ... 19
表3.2 模式建立所採用各颱洪事件之資料起迄時間 ... 22
表3.3 不同臨前序列長度下聖帕颱風預報一小時與觀測值之間平均誤差 ... 23
表3.4 不同臨前序列長度下柯羅莎颱風預報一小時與觀測值之間平均誤差 ... 23
表4.1 基隆河流域各集水分區面積、河川長度、平均坡度 ... 26
表4.2 雨量站座標位置 ... 28
表4.3 率定之係數值 ... 29
表4.4 率定之係數值 ... 30
表4.5 搜尋最佳參數 ... 32
表4.6 大華橋洪水位各颱風預報一小時誤差表 ... 34
表4.7 最小二乘法與結合 SA 比較平均誤差(M) ... 35
表4.8 聖帕颱風於基隆河各水位站預報一小時誤差表 ... 38
表4.9 柯羅莎颱風於基隆河各水位站預報一小時誤差表 ... 41
表4.10 水壓式水位計 ... 45
表4.11 雷達波式水位計 ... 45
表4.12 員山子觀測與計算分洪量之誤差 ... 47
表4.13 員山子分洪日期、事件名稱及分洪量 ... 47
表A-1 河川不恆定流模式採用河段斷面說明(新店溪、大漢溪部份) ... 73
表A-2 河川不恆定流模式採用河段斷面說明(基隆河部份) ... 74
圖目錄
圖2.1 模擬退火演算法物理意義圖解(童慶斌,2008) ... 6
圖2.2 模擬退火演算法搜尋過程圖解(許家源,2003) ... 7
圖2.3 模擬退火法之搜尋求解示意圖 ... 9
圖2.4 模擬退火演算法步驟流程圖 ... 11
圖2.5 移步機率函數圖形 ... 12
圖3.1 球與直角座標 ... 20
圖3.2 集水區參數求鄰近解流程圖 ... 20
圖3.3 模式建置流程圖 ... 24
圖4.1 基隆河流域圖 ... 25
圖4.2 基隆河河道分配圖 ... 27
圖4.3 基隆河大華橋集水區徐昇多邊形網 ... 28
圖4.4 大華橋參數率定 ... 31
圖4.5 辛樂克颱洪事件大華橋水位站 ... 31
圖4.6 聖帕颱洪事件大華橋水位站 ... 32
圖4.7 柯羅莎颱洪事件大華橋水位站 ... 33
圖4.8 薔蜜颱洪事件大華橋水位站 ... 33
圖4.9 大華橋洪水位各颱風預報一小時誤差 ... 34
圖4.10 最小二乘法與結合 SA 比較平均誤差圖 ... 35
圖4.11 聖帕颱風五堵水位站 ... 37
圖4.12 聖帕颱風江北橋水位站 ... 37
圖4.13 聖帕颱風社后橋水位站 ... 38
圖4.14 聖帕颱風基隆河各水位站洪水位預報一小時誤差 ... 39
圖4.15 柯羅莎颱風五堵水位站比較圖 ... 39
圖4.16 柯羅莎颱風江北橋水位站比較圖 ... 40
圖4.17 柯羅莎颱風社后橋水位站比較圖 ... 40
圖4.18 柯羅莎颱風南湖大橋水位站比較圖 ... 41
圖4.19 柯羅莎颱風基隆河各水位站洪水位預報一小時誤差 ... 42
圖4.20 員山子分洪設施位置 ... 43
圖4.21 員山子分洪設施名稱 ... 44
圖4.22 員山子分洪設施示意圖 ... 45
圖4.23 柯羅莎颱洪事件員山子分洪量 ... 46
圖4.24 薔蜜颱洪事件員山子分洪量 ... 46
圖4.25 柯羅莎颱洪推估上游入流量 ... 50
圖4.26 員山子攔河堰前後水位示意圖 ... 50
圖4.27 柯羅莎颱洪攔河堰前後水位 ... 51
圖4.28 柯羅莎颱洪員山子各測站水位 ... 52
圖4.29 柯羅莎颱洪分洪堰前水位 ... 53
圖4.30 柯羅莎颱洪計算分洪流量 ... 53
圖A-1 互相交叉之特徵線示意圖 ... 62
圖A-2 第二種多方式特徵法(MMOC-II) ... 62
圖A-3 程式架構及資料連結展示 ... 66
圖A-3 程式架構及資料連結展示(續) ... 67
圖A-4 本模組 CCCMPORT.A01 輸出示意圖 ... 71
圖A-5 本模組 CCCMPORT.A06 輸出示意圖 ... 71
圖A-6 本模組 CCCMPORT.A07 輸出示意圖 ... 72
第一章 緒論
1.1 前言
台灣全島雨量豐沛,年平均降雨量為世界平均值的2.5 倍,高達 2500 公釐。
降雨主要是集中在夏季,其每年平均約三至四個颱風,所帶來的豪雨占全年降雨 量的大多數。由於颱風夾帶的雨量大且速度快,若不能立即宣洩排洪,則勢必造 成災害。
台灣的主要河川大部分都發源自中央山脈,分別向東西方注入太平洋或台灣 海峽,由於中央山脈位置及地勢偏東的緣故,故台灣的河流東短西長。全國大小 河川共計有 129 條,長度均甚短,坡陡水急。淡水河為台灣第三大河川,其支流 新店溪、基隆河、大漢溪等河流的氾濫平原,造就今日台北盆地,就像世界上的 古文明發祥地黃河、中東的兩河流域一樣,都是起源於氾濫平原上。因為那裏土 壤肥沃,灌溉、飲水所需水源取之方便。淡水河流域包括台北市、台北縣大部份、
桃園及新竹兩縣之小部份;其中台北市二百多萬人,工商業發達,為我國政治及 經濟文化中心。
綜合以上所述,台灣河川均短且急,地形大部分屬於高山與丘陵地。由於降 雨季節集中,雨量豐沛,且經過都市化的開發後,不但地表逕流量增加,更讓集 流時間大幅縮短,導致洪水來不及宣洩,往往在短時間內使得河川水位驟升,造 成溢堤,淹水災害頻傳,附近居民苦不堪言。故如何及時且有效率的預測降雨量 與水位的關係,建立具有洪水預警系統的模式,能夠在災害發生前進行疏散驅離,
將損失降至最低,此乃重要研究方向之一。
1.2 文獻回顧
1.2.1 降雨逕流或水位推估
關於降雨-逕流推估方法,最早可以回溯自 Sherman 在 1932 年所發展的單位歷 線之理論,此理論建立於假設符合線性系統,即在某特定單位降雨延時內可產生 一條逕流歷線,在不同強度的降雨但相同延時內所產生之逕流歷線,可以該單位 歷線正比關係進行推估。因此,水文工程師可利用集水區過去的水文歷史記錄資 料計算其單位歷線,應用假設線性的特性,以疊加與正比關係,推算集水區現在 發生暴雨所產生的逕流歷線。Nash 在 1957 年提出將集水區視為由多個串聯的線性 水庫所組成,假設每個水庫之出流量與其蓄水量呈線性正比,而最後一個水庫出 流量即為此集水區之瞬時單位歷線。Henderson and Wooding 在 1964 年應用運動波 理論,推導集水區漫地流集流時間,配合Wooding 於 1965 年提出的渠流部分之集 流時間,建立一套漫地流及河川逕流量演算方法,藉以模擬降雨和集水區出口逕 流量之關係。日本科學技術廳菅原正已(Sugawara)於 1971 年所推廣之水筒模式法,
其概念乃是將集水區逕流機制,代換成用數個貯留型之模型容器所組合,考慮集 水區內複雜之水文因子,如入滲、蒸發、貯留、地表逕流、中間流、基流及滲漏 損失等現象,推導出降雨和逕流量之關係。近年來由於電子計算機計算能力之提 昇,類神經網路廣泛的應用在推估降雨-逕流上,其相關論文皆陸陸續續地發表出 來,Zhu and Fujita 於 1994 年應用模糊推論的方式,採用簡單三角形,定義降雨及 逕流的隸屬函數及建立降雨-逕流之模糊關係,以預測下一刻之逕流量。Komda and Makarand 於 2000 年以倒傳遞類神經網路進行洪水位之預測。Lin and Chen 於 2005 年藉由應用全面監督式訓練法,建立輻狀基底函數之網路架構,以作為預報洪水 流量之模式。吳南靖等於2004 年應用序率方法推估河川特定點進行水位預報之研 究。
1.2.2 模擬退火演算法
最早發展模擬退火( Simulated Annealing,SA)之優選模式的是人們稱之為蒙地 卡羅(Monte Carlo)演算法,是 Metropolis 等人在 1953 年提出,不過當時處理的 是 提 供 高 速 電 腦 處 理 有 限 空 間 最 大 放 置 問 題 (two-dimensional rigid-sphere system)。直到Kirkpatrick 等人在 1983 年提出並成功地應用在大尺度組合最佳化問
題中,它是蒙地卡羅演算法的推廣。此種隨機搜尋之演算法,將求解流程比擬為 熱處理中退火之現象,退火是一種物理過程,當金屬物體加熱至一定的溫度後,
它的所有粒子在狀態空間中具高活性並呈現自由運動。隨著溫度的下降,這些粒 子逐漸停留在不同的狀態下。當溫度降至最低時,粒子重新並且以最佳的結構排 列,而粒子也就是以波茲曼(Boltzamnn)概率進行分布,因此如何設定溫度作有 效率的搜尋就變成整個模擬退火法最重要的一環。其相關研究:Kirkpatrick(1983)
應用於大尺度組合優選之問題上;Huang and Romeo(1986)介紹其降溫梯度參數 之設定;Dougherty and Marryott(1991)介紹最大搜尋次數與決策變數的關係;
Mayer et al.
(1998)農場管理問題
;童慶斌等(2000)應用模擬退火法在計算集水 區平均雨量之應用等等。雖然此演算法已經有五十幾年的歷史,但基於它的簡單 性和實用性,所以仍被廣泛地運用在各種最佳化問題中。1.3 研究目的
基隆河為大台北地區民眾的生活經濟動脈,但因都市化與工商業急速的發 展,人口密度不斷的增加,導致近年來基隆河兩岸及其流域內土地過度被開發利 用,降雨逕流亦隨之增加。
基隆河於關渡附近匯流於淡水河,其中、下游河道蜿蜒曲折而平緩。由於部 分河段之河寬過於狹窄,因此基隆河本身自然排洪條件就不佳,每逢發生豪大雨 時即傳出淹水災害。尤其是 1987 年 10 月琳恩颱風所帶來的豪雨,造成台北市部 分行政區及台北縣汐止、瑞芳鎮及基隆市等廣大地區遭受水患,土地淹沒,損失 慘重。近年來,政府致力於進行各項工程建設,台北地區防洪計畫大致完工使用,
但淹水事件卻仍層出不窮,如 1996 年賀伯颱風及 2000 年象神颱風等等,皆造成 瑞芳、汐止及南港等地區發生嚴重淹水災情。有鑒於此,行政院指示:「宜考量優 先推動員山子分洪工程之可行性。」,並於2000 年 11 月行政院院長裁示立即推動 員山子分洪計畫之相關工作。
為使基隆河有效減低中、下游地區之洪水負擔,於2002 年台北縣瑞芳鎮境內 興建員山子分洪道設施,將颱風或豪大雨造成的洪水進行部分疏導,而原有的水 位預報模式(吳等,2004)計算大華橋水位站,便產生高估水位的現象。本研究目 的,在於考量分洪道之分洪流量並且納入每年洪水沖刷或淤積河道之影響等因素
修正模式,使大華橋水位能夠更精準的預報。
本模式應用柯羅莎(Krosa)、薔蜜(Jangmi)等颱風預報大華橋水位之變化,皆 獲得良好之驗證。其預報結果可提供全流域不恆定流河川水位預報模式之基隆河 上游水位邊界點條件,對於下游洪水預警,作更加準確的預報,以爭取沿岸居民 疏散時間。
第二章 模擬退火演算法之理論
近年來由於「最佳化」一名詞在各種不同領域中快速地廣泛使用,例如電機 工程、電腦科學、工業工程以及物理研究等。如何在茫茫未知解的空間中,尋找 出最佳的解,便成為重要課題之一。
面對這一類複雜優選相關問題時,一些啟發式的方法便開始逐漸發展成型,
目的就是希望利用電腦作較有效率或效能的搜尋,在合理的時間內找出最佳解。
通常在求解最佳化一系列的步驟中,即為一種演算法。當面臨需要使用演算法求 解問題時,則須將所遇到之問題模式化,設定邊界條件及目標函數式可將問題公 式化,為因應各種不同的問題與需求,其所設定的目標函數與邊界條件也都有所 不同,考量問題之規模與尺度,使用最合適的演算法,將可以最有效率地達到最 佳化之目的。
啟發式演算法主要是利用人工智慧技術,搜尋結構複雜問題之最佳解,經常 可以有效率地求得可被接受的答案。禁忌搜尋法(Tabu Search) 採用傳統搜尋技術 並且結合記憶之功能,模擬退火法(Simulated Annealing Algorithm)以模仿晶體高溫 下逐漸冷卻至達到穩定解之過程,遺傳演算法(Genetic Algorithm)模仿生物適者生 存之法則,透過基因有擇優、交配及突變的能力繁衍後代以適應周遭環境以及蟻 行演算法(Ant Algorithm)模仿螞蟻行進以最短路徑之特性搜尋最佳解。上述這些近 數十年發展的啟發式搜尋方法,不外乎為聰明的高等試誤法,其利用模擬大自然 之各種特性,設計不同的搜尋規則,建構出更加具有求解最佳解能力的演算方法,
以因應未來高科技發展可能帶來更加複雜的最佳化問題。(上述之內容, 主要參考 自童慶斌,2008)
2.1 模擬退火法之原理及物理現象
根據物理概念,在加熱金屬固體物質時,隨著溫度不斷提高,粒子
將會因為能
量高、活性大,而使得
固體內部粒子運動愈趨激烈,其粒子分布與原始位置不斷偏離。當溫度升高至超出固體熔點溫度時,其晶體結構之規則性將被破壞,熔化成 液態。粒子彼此間之結構排列從有順序規則之固態轉變為沒有順序之液態,這種 過程稱之為熔解。熔解過程主要目的在於消除固體中原本可能存在之不均勻狀 態,高溫使粒子有機會重新排列形成各種不同的型態,隨後進行冷卻過程使其運 動逐漸趨緩,也減少了粒子可能型態的數目。溫度逐漸降低,粒子之運動愈來愈 有規則。當溫度降至結晶溫度時,系統最終被冷卻到基底狀態,粒子之運動則會 圍繞著晶格作微小震動,代表粒子之間能量與晶格排列達到最佳狀態。液體凝固 成為固態之晶體,稱之為退火。所以系統的最高溫度要夠高且降溫過程不能夠太 慢,才算是物體完整地退火過程(李世炳、鄒忠毅,2002),其上述之退火結晶過程 如圖2.1。
圖2.1 模擬退火演算法物理意義圖解(童慶斌,2008)
模擬退火法便是模擬液體的退火現象,相對應於物理系統的退火過程,解的 集合就好像是物質的型態,不同的結晶狀態代表著不同的可行解,當液體重新結 晶成固體最穩定的型態時,也就是搜尋至最佳解。在高溫狀態下物質內的粒子擁 有高活性,可任意重組液體型態,而演算法之原始溫度設定於高溫狀態便類似此 種情況,因其最佳化之參數擁有高度的可變性,目的為使系統之初期目標函數較 不受目標函數值優劣影響,可有各式各樣不同的組合型式,表示搜尋的過程中其 解可能跳出最佳解區域或跳入較劣解之區域。在緩慢降低溫度的過程中,粒子活 性隨時間漸漸減弱,熔融狀之液體將慢慢結晶成為穩定固體型態,直到最後溫度 冷卻至固態穩定晶體析出,而演算法溫度慢慢下降時,解集合也漸漸無法接受任 意組合之型式,而開始傾向往表現出目標函數值較優良的方向搜尋,接受其解集 合,直到系統溫度達到設定之最終溫度則停止搜尋,其解集合也將會收斂至一極 佳值(童慶斌,2008),其搜尋系統如下圖 2.2 所示。
圖2.2 模擬退火演算法搜尋過程圖解(許家源,2003)
不過與物理退火現象較為不同之處為液體從高溫熔融退火至低溫結晶體析 出,其最後之結晶型態是為最穩定狀態,而在最佳解之搜尋系統中,從系統高溫
狀態下搜尋解集合,隨著系統設定之降溫,其搜尋解最後收斂於一極佳值,但此 極佳值有可能是為區域最佳解(Local Optimal)而非系統中之全域最佳解(Global Optimal)。利用模擬物質退火之過程中粒子不斷重新移動而改變排列組合,如同求 解最佳化問題之演算程序,因此整理如表2.1 顯示二者之間進行類比對照關係表。
表2.1 物理系統與搜尋最佳解類比對照表
物理系統 搜尋最佳解
型態(State) 可行解(Feasible Solution)
內部能量(Energy) 目標函數(Objective Function)
粒子活性(Particle Activity) 參數可變性(Parameter Variability)
退火(Annealing) 模擬退火(Simulated Annealing)
基態(Ground State) 最佳解(Optimal)
2.2 模擬退火法之搜尋方式及演算步驟
模擬退火法屬於區域搜尋法(Local Search Method)之ㄧ,結合了最陡坡降法及 隨機搜尋的方式來尋找目標函數之最佳解,其移動規則為鄰近空間中隨機挑選一 個較佳之鄰近解,並將此相對較佳之鄰近解取代為目前解。模擬退火法計算開始 時即在可行解之集合中隨機產生一個初始解,利用鄰近區域搜尋較佳之結果取代 目前較劣之結果,但容易使求出之解陷入區域最佳解,所以模擬退火法在區域搜 尋之過程中,透過機率數值決定是否可以往較高能量方向移動,而此決定機制即 為波茲曼機率分佈(Boltzmann probability distribution)。運用機率方式允許接受較差 之劣解,藉此跳脫求解過程中陷入區域最佳解之困境,進而逐步逼近至全域最佳 解(張,2003),圖 2.3 為顯示模擬退火法之搜尋求解示意圖。
圖2.3 模擬退火法之搜尋求解示意圖
其相關名詞之涵義如下:
(1) 可行解:所謂『解』,是代表一種狀態,在模擬退火法則為每個決策變數 值,其滿足問題中所有限制之條件稱為可行解。
(2) 初始解:從可行解空間中隨機產生初始點作為起點,此解即為初始解。
(3) 目前解:搜尋過程之目前所在點。
(4) 鄰近解:目前解之鄰近區域中任一可行解,此鄰近區域可由演算法依問 題特性而定義。
(5) 區域最佳解:鄰近區域中最佳解。
(5) 全域最佳解:整個可行解區域中最佳解。
模擬退火法包含了 Meteroplis 演算法及退火程式兩部份。基本內部迴圈乃是 以Meteroplis 演算法為主,降溫機制由退火程式作為外部迴圈,即在每次溫度降低 時執行Meteroplis 演算法,其搜尋最小值演算步驟分別敘述如下,圖 2.4 所示為其 演算步驟流程圖。
(1) 設定初始狀態,如在行解空間中隨機選定初始解X ,做為現行解0 Xi = X0,並 計算其目標函數值F X( )i =F X( 0)、初始溫度Ti = 、馬可夫鏈長度(T0
每一溫度
階段的搜尋次數
)等。(2) 依照目前解X ,隨機選取鄰近解i Xi+1,並計算其目標函數值F X( i+1)。 (3) 計算目標函數值差值,Δ =F F X( i+1)−F X( )i 。
(4) 判斷是否接受移步:
1 0
P( )
exp( ) 0
i
if F
F F
T if F
⎧ Δ <
Δ = ⎨ ⎪ ⎪ −Δ Δ ≥
⎩
○1 若
Δ < F 0
,表示鄰近解Xi+1優於目前解X ,則接受此次鄰近解為下次目前解。 i ○2 若Δ ≥ F 0
,表示鄰近解Xi+1比目前解X 差,用波茲曼機率分佈iP exp( )
i
F T
= −Δ
,並產生一隨機數值(介於 0~1)R,判別是否接受此鄰近解Xi+1,若P ≥ R 則接 受本次鄰近解Xi+1為下次目前解X ,反之則否,依然以本次目前解i X 當作下次 i 移步基準之目前解。
(5) 判斷條止條件(收斂條件),若符合則停止演算法,並以目前解作為最佳解,若 否,則進行步驟(6)。
(6) 降低溫度,Ti+1= × 產生下一個溫度,其中k Ti
k
介於0~1 之間,重複步驟(2)繼續搜尋。
圖2.4 模擬退火演算法步驟流程圖
圖2.5 移步機率函數圖形
2.3 模擬退火演算法影響因子之探討
在模擬退火法中,各項參數值建立對最佳解之搜尋皆會直接影響其品質,若 是發生設定不恰當,即有可能造成搜尋時間上的增加或是降低最佳解之品質,因 此以下將針對各項參數值設定及影響因數作討論。然而並非同樣的退火時間表都 適用於每一種問題,須視其問題之尺度與規模而決定退火參數為何,以下為前人 針對退火參數設定之建議,但其適合之退火時間表,仍需要依照問題重新定義。
1. 初始溫度設定
起始溫度設定在模擬退火法中扮演著重要的角色,因溫度與目標函數差值變 化,決定劣解被接受機率之大小,因此溫度的變化在退火演算法中,佔有相當重 要的影響力,本研究利用波茲曼機率分佈公式(Boltzmann distribution),來決定是否 接受劣解。在整個搜尋過程中,從波茲曼機率分佈
P exp( )
i
F T
= −Δ
可得出,機率基本上與溫度T 及目標函數之差值i Δ 相關,其中固定溫度F T 的條件下,劣解被接受的i 機率為P,與目標函數差值Δ 呈反比關係。一般來說,在退火程式開始前部分,F 起始溫度T 必需要夠大以確保搜尋初期時,對於所有鄰近解,有著相同的挑選機0 率,其溫度較高被接受機率值也相對提高,允許接受劣解機率較大,以避免落入 區域最佳解,雖然過高的溫度設定會造成計算上的浪費,不過設定高溫度的結果 普遍來說都優於設定過低的結果(Dougherty and Marryott, 1991),所以一般設定適
合溫度,就是讓初期鄰近解被接受率大約在 80%以上,期望搜尋過程中不要太快 收斂於某一區域最佳解中。隨著溫度降低至接近凝固點時,理論上此時已逼近全 域最佳解,達到最佳狀態,進而機率也隨著溫度而降低,無須再逆向搜尋。所以 在以模擬退火法解決一最佳化求解之問題時,不旦要定義出適當的限制式與目標 函數式,更要適時地調整其退火參數以及最重要之核心方程式,才能有效率地解 決問題(童慶斌,2008)。
2. 馬可夫鏈
馬可夫鏈的長度(length of Markov chain)(L)即為每一個溫度階段的最大搜尋次 數。模擬退火法可視為多階層之馬可夫鏈,每一溫度階段的搜尋次數長短,也將 影響到最佳解之品質,Kirkpatrcik et al. (1983) 建議搜尋次數為決策變數個數的倍 數;Arts and van Laarhoven (1985) 曾經證明只要搜尋次數夠多,模擬退火法最後 找到的最佳解,為全區域最佳解的機率可以高達百分之百,除此之外,一些前人 所建議每一溫度階段之搜尋次數:Huang et al. (1986) 則採用隨著不同溫度而改變 搜尋次數,當溫度愈降愈低時,接受鄰近解的機率已經大大的不如初始的接受機 率,便不再需要初期溫度階段的搜尋次數,而可以隨著溫度的降低減少搜尋次數,
以免造成效率上的浪費;Dougherty and Marryott (1991),實際將搜尋次數訂為決策 變數個數的 100 倍。不過,至今相關的報告,大多無法肯定共通的準則來決定每 一溫度階段的搜尋次數,還是要視問題的尺度大小與特性而定。
3. 降溫梯度
降溫梯度
k
是退火時間表中控制溫度參數從初始高溫度經過退火降溫之過 程,最後到達最終溫度Ti → 的參數(童慶斌,2008),其關係式如式(2.1)所示: 0Ti+1 = ×k Ti , 0< <k 1 (2-1) 其中Ti+1代表降溫後溫度,T 代表降溫前溫度值,i
k
代表每次溫度改變速率。上式2.1 中,其降溫梯度(
k
)為 0~1 之間係數值,而降溫梯度的高低會影響到 求解品質與其運算時間;由前一節中,可瞭解 Metropolis 程式接受鄰近解的機率 程度受到溫度影響;k
值越大時,代表降溫速度越趨於緩慢,於整個退火過程有較 高接受率,能有足夠的空間搜尋,反之則代表降溫速度越快速,接受率較低,容 易導致陷入區域最佳解中。在各研究中,皆曾對於溫度下降速率與機制進行探討。一般而言,設定
k
值介於0.80~0.99 之間 (Kirkpatrick et al.,1984),整個退火過程 皆有不錯的搜尋表現。不過降溫梯度與每一溫度階段的搜尋次數是有相關的,當搜尋次數設定不夠 時,降溫梯度則會變的更加重要,也就是說可以將降溫梯度設定為較大的值(趨近 於 1),讓每一次降溫幅度不至於太大,不然的話,太快的降溫再加上每一個溫度 階段中搜尋次數又不夠多,這樣子所搜尋解空間範圍有限,無法有效的搜尋,而 且最後所收斂的結果也有可能不夠理想,所以在電腦運算合理時間內,設定最佳 的降溫梯度加上合適的馬可夫鏈長度,才能夠有效率的搜尋出最佳解(童慶斌,
2008)。
4. 停止條件
停止條件是判斷是否該終止演算法之搜尋,若設計不適宜,可能導致演算過 程中,結束時間過早或執行時間過長。一般設定的停止條件方式如下(Dougherty and Marryott, 1991):
(1) 達到終止溫度,當溫度低於終止溫度則停止演算。
(2) 當連續降溫 K 次始終無法改善最佳解,則停止演算。
(3) 當搜尋解之數量已達到最大迭代次數,則停止演算。
第三章 模式建立
基於河川水位具有連續性的變化,即現時刻水位與前一段時刻水位記錄有極 大之關聯性,所以採用序率方法來推算水位變化乃可行之道。然而,影響水位之 因子,除了臨前水位變化之外,尚有其集水區內的降雨量、員山子分洪工程等自 然或人為因素。故本研究以序率方法為基礎 ,藉由最小平方法得出一條具遞迴關 係特性的預報方程式,最後應用模擬退火法修正方程式之參數,使其更佳適應逐 年河道因洪水而造成的改變,並達到準確預報之目的。
3.1 序率分析法
假設河川水位站在某個時刻的水位,為其水位時間序列記錄以及若干個影響 變因之函數。其影響變因係為集水區雨量站及員山子分洪水位亦為時序列之資 料,且考量影響變因對該處之水位變化可能有立即影響,其關係寫成函數式表示 如下:
( , , )
L= f L
ξ
h (3-1)其中L 、
ξ
與h
各為一數學向量,表示如下:{ }
{ }
{ }
( ), 1, 2, ( ), 0,1, ( ), 0,1, L L t i t i
t j t j h h t k t k
ξ ξ
= − Δ =
= − Δ =
= − Δ =
(3-2)
Δ t
為水位或流量觀測記錄之時間間距,本研究設定時間間距Δ t
為一小時,L 代 表該水位站臨前狀態,h
為員山子水位高,而ξ
為 M 個項目所組成的矩陣,表示 如下:{
1 2}
( ) t ( ), ( ), , t t
M( ) t
ξ = ξ ξ ξ
(3-3)ξ
之項目,其變因主要為集水區內或鄰近雨量站當時及臨前狀態之降雨量。在鄰近雨量站方面,吾人擬先採用基隆河流域的雨量站,利用徐昇氏多邊形 法將集水區分區,若以大華橋水位站上游之徐昇氏多邊形有包含到集水區,則該 雨量站就選入,否則就不列入考慮。
本模式所採用之降雨量,有別於降雨-逕流模式所採用扣除截流、窪蓄、入滲 之後所得的有效降雨,而是雨量站得直接降雨記錄或預報降雨。同時,本模式無 須另行計算各子集水區在各時間步的平均有效降雨,免除推估有效降雨及集水區 平均雨量所可能遭遇的不確定因素。
除此之外,降雨在空間上的變異性對河川水位變化之影響也能夠有所反應。
實際上,本模試考量臨前時序列水位狀態已隱含基流量、截流量、窪蓄量、入滲 量等因素對河川水位之影響。
由於各因素對河川水位之影響隨時間越遠而越小,吾人可進一步假設一個影 響時間長度,超過這個長度其影響性可忽略,則其函數關係可表示成:
{ }
{ }
{ }
( ), 1, 2, , ( ), 0,1, , ( ), 0,1, , L L t i t i N
t j t j N h h t k t k N
ξ ξ
= − Δ =
= − Δ =
= − Δ =
(3-4)
其中 N 為一個待決定的正整數,表示臨前記錄之影響時間長度。本研究設定 N 為 6,換句話說,也就是臨前記錄之影響時間為六小時,其不同臨前記錄影響長 度如表3.3、3.4 所展示。
假設(3-4)式的函數關係為各變因之線性組合,可表示成:
0
1 1 0 0
( )
N i( )
M N ij i( )
N k( )
i i j k
L t a a L t i t b ξ t j t c h t k t
= = = =
= + ∑ − Δ + ∑∑ − Δ + ∑ − Δ
(3-5) 其中,a 、0 a 、i b 、ij c 為待率定的係數。 k這些係數可從歷史洪水記錄中透過最小平方法率定而得。此函數關係之方程 式雖為線性表示式,但因具時序性,變因與變量間具遞迴關係,可表現非線性之 物理現象,且具有隨資料不斷更新而即時修正之功能。
本模式在應用時,當掌握到所有影響河川水位的因素(例如:降雨量、員山子 水位甚至上游水庫…等)時,透過水位的時序遞迴關係,即可預報未來一段時間的 水位變化。
3.2 最小二乘法 式(3-5)可改寫成
0 1 1 2 2 z z
L X= +X A +X A + +X A (3-6)
其中Z = +N (N+ ×1) (M + ,上式亦可改寫成: 1)
0 1 11 2 21 1
1
0 1 12 2 22 2
2
0 1 1 2 2
z z
z z
n n n z zn
X X A X A X A L
X X A X A X A L
L X X A X A X A
+ + + +
⎡ ⎤
⎡ ⎤ ⎢ ⎥
⎢ ⎥ ⎢ + + + + ⎥
⎢ ⎥ =⎢ ⎥
⎢ ⎥ ⎢ ⎥
⎢ ⎥ ⎢ ⎥
⎢ ⎥ + + + +
⎣ ⎦ ⎣ ⎦
(3-7)
如以矩陣來表示,則可改寫成:
0
1 11 21 1
1
2 12 22 2
2
1 2
1 1 1
z z
n n zn
n
z
L A A A X
L A A A X
X
A A A
L X
⎡ ⎤
⎡ ⎤ ⎡ ⎤ ⎢ ⎥
⎢ ⎥ ⎢ ⎥ ⎢ ⎥
⎢ ⎥ ⎢= ⎥ ⎢ ⎥
⎢ ⎥ ⎢ ⎥ ⎢ ⎥
⎢ ⎥ ⎢ ⎥ ⎢ ⎥
⎢ ⎥ ⎣ ⎦
⎣ ⎦ ⎢ ⎥
⎣ ⎦
(3-8)
[ ]
nL
×1 =[
n z× +A
( 1)] [
(zX
+ ×1) 1]
可簡寫成 L AX= 兩邊各乘A : T
T T
A L A AX= (3-9)
兩邊再各乘以(A A 的反矩陣,即 T )
1 1
(A AT )− A LT =(A AT ) (− A A XT ) (3-10)
故得X =(A AT )−1A LT
3.3 模擬退火演算法於集水區參數調整之運用
1.變數個數之選定
在模擬退火法中,各項參數值建立對最佳解之搜尋皆會直接影響其效率,若 是發生設定不恰當,極有可能造成搜尋時間上的增加或是降低最佳解之品質,因 此參數的選定,實乃重要。除此之外,選定變數的多寡也是影響效率的一大因素,
應用吾人在上一節中藉由最小平方法所求出的X 向量,從中挑選出數個變數值,
但若挑選得過多,則會造成電腦搜尋時間上的浪費;若挑選得過少,則會造成行 解空間上搜尋不足,導致最佳解品質降低,並且挑選出得變數值須具有代表性,
所謂的代表性,在此是指經由修正挑選之變數,可以針對每年洪水所造成河道之 變化依然進行準確預報的功能。
因此本文設定挑選三個變數值,分別為X( )
α
、X( )β
、X( )γ
,其中α 、β
、γ
為代表員山子水位站的t − Δ 2 t
、t − Δ t
、t 時刻。選擇員山子水位站作為變數,主 要是考慮到員山子水位站位於基隆河的主河道上游處,能有效的反映給大華橋水 位站並且當河道上游洪水流量足以啟動員山子分洪機制時,表示其河川水位高於 標高63 公尺之分洪堰,而龐大的流量則極有可能造成下游河床的改變,故挑選員 山子水位站參數作為變數值。2.目標函數
本優選模式之目標函數為令計算水位與觀測水位差值之平均誤差為最小,該 目標函數的訂定,主要目的是使計算水位歷線儘可能接近觀測水位歷線。可表示 如下式(3-7):
1
Z n ( ( ) ( ))
i
Min A i X L i
=
=
∑
− (3-11)其中A i : A 矩陣的第 i 列,為影響變因的時間序列。 ( ) ( )L i : L 向量的第 i 值,為觀測水位。
X :藉由最小平方法所求出的一數學向量。
3.設定相關參數
在模擬退火法中,欲搜尋之變數通常都是隨機亂數給定其初始值,但在本研 究中變數初始值是藉由最小平方法計算所得出的值當作其初始值。而所需設定之 其它相關參數包括初始溫度、馬可夫鏈長度、降溫梯度、最終溫度、停止條件分 別如下表3.1 所示:
表3.1 相關參數設定
參數 設定值
初始溫度 10℃
馬可夫鏈長度 10000
降溫梯度 0.1
最終溫度 10−10 C
停止條件 階段溫度低於最終溫度時停止搜尋
4.搜尋方式
本研究之搜尋方式,採用類似於球型之搜尋法,如圖 3.1(左)所示,有別於傳 統方式一次只能搜尋單一方向,其優點為變數一次可以全部變動,適合於整組變 數作調整。首先,先決定球半徑r 為5 10× −9的移步量,再將θ 與
φ
各隨機產生八組 角度,然後隨機挑選其中一組計算其X( )α
、X( )β
、X( )γ
,作為鄰近解。式(3-12) 為其直角座標方程式。( ) sin cos ( ) sin sin ( ) cos
X r
X r
X r
α θ φ
β θ φ
γ θ
⎧ =
⎪ =
⎨ ⎪ =
⎩
(0≤ ≤
θ π
and 0< ≤φ
2 )π
(3-12)圖3.1 球與直角座標
上述之整個集水區參數應用於模擬退火演算法鄰近解的流程圖,如圖 3.2 所 示:
start
圖3.2 集水區參數求鄰近解流程圖
3.4 建置步驟
若假設L 函數為各變因之線性組合,表示為:
0
1 1 0 0
( ) ( ) ( ) ( )
N M N N
i ij i k
i i j k
L t a a L t i t b
ξt j t c h t k t
= = = =
= + ∑ − Δ + ∑∑ − Δ + ∑ − Δ
(3-5)利用線性迴歸之方法,以最小平方法求出各個係數(各係數矩陣的各元素),則 達到建立模式之目的。模式建立之步驟如下:
(1)選定基隆河流域集水區內各個雨量站。
(2)篩選雨量站:
○1 在基隆河集水區內繪製徐昇多邊形網。
○2 在各個徐昇多邊形網中,若有包含到上游集水區河川水位站者,其該雨量站 之降雨量則判定為影響河川上游水位之雨量因素。
(3)判斷是否還有其他影響河川水位之因素,如:員山子水位高、上游水庫…等。
(4)將歷史颱洪事件中的觀測記錄,分成率定組及驗證組兩部分統計樣本,以茲比 較,如表 3.2 所示。
(5)暫定影響河川水位之臨前序列長度,應用最小平方法,從率定組統計樣本資料 求出預報方程式中的各係數。
(6)從驗證組中,挑選其兩場颱洪事件,作為各觀測資料與預報方程式輸出值作誤 差方均根之計算。
(7)嘗試不同的影響河川水位之臨前長度,進行步驟(5)及(6),如表 3.3 所示,挑選 出最佳的臨前序列影響長度。
(8)重新推求參數:利用前一步驟所求得的臨前序列影響長度,重新推求預報方程 式中的各係數。
(9)從各係數中挑選員山子水位之三個係數值,分別為X( )
α
、X( )β
、X( )γ
, 其 中α 、β
、γ
為代表其員山子水位的t − Δ 2 t
、t − Δ t
、t 時刻。(10)從率定組統計樣本中,應用模擬退火法求得最佳參數值。
(11)模式使用新的一組參數,進行計算驗證組。
表3.2 為模式建立所採用各颱洪事件之資料起迄時間表,分成率定組與驗證組 兩組。率定組之所以蒐集多場颱洪歷史事件記錄,主要是針對員山子分洪流量,
只有當基隆河流域上游地區其降雨量大到足以啟動員山子分洪機制時,才會有分 洪情況產生,由於每場颱洪事件其分洪記錄大部分都屬於一小段時間,所以必須 採用多場颱洪事件來彌補分洪流量記錄不足之地方。
表3.2 模式建立所採用各颱洪事件之資料起迄時間
組別 颱洪事件 時間
率定組
瑪莎 2005/08/04 20:00 ~ 2005/08/05 20:00 泰利 2005/08/31 22:00 ~ 2005/09/02 06:00 龍王 2005/10/02 06:00 ~ 2005/10/02 12:00 豪雨 2006/09/10 16:00 ~ 2006/09/11 03:00 帕布 2007/08/07 22:00 ~ 2007/08/08 23:00 韋帕 2007/09/18 08:00 ~ 2007/09/19 01:00 米塔 2007/11/26 09:00 ~ 2007/11/27 23:00 辛樂克 2008/09/12 01:00 ~ 2008/09/16 00:00
驗證組
聖帕 2007/08/17 07:00 ~ 2007/08/21 00:00 柯羅莎 2007/10/06 06:00 ~ 2007/10/07 05:00 薔蜜 2008/09/28 06:00 ~ 2008/09/29 12:00
表3.3 不同臨前序列長度下聖帕颱風預報一小時與觀測值之間平均誤差
N MAE N MAE
1 0.088 6 0.050
2 0.054 7 0.184
3 0.052 8 0.154
4 0.054 9 0.156
5
0.049
10 0.154上表中MAE:
1
1 n ˆi i
i
z z n =
∑
− ,其中 ˆz 為計算值;i z 為觀測值;in
為計算樣本數表3.4 不同臨前序列長度下柯羅莎颱風預報一小時與觀測值之間平均誤差
N MAE N MAE
1 0.161 6 0.107
2 0.115 7 0.173
3 0.095 8 0.135
4 0.096 9 0.113
5
0.090
10 0.114上表中MAE:
1
1 n ˆi i
i
z z
n
∑
= − ,其中 ˆz 為計算值;i z 為觀測值;in
為計算樣本數β γ
圖3.3 模式建置流程圖
第四章 模式驗證
4.1 研究區域
基隆河發源於台北縣平溪鄉菁桐山,先後由鰈魚坑溪、東勢坑溪、大武崙溪、
瑪陵坑溪、鹿寮溪、保長坑溪、北港溪、大坑溪等大小支流匯合,流經平溪鄉、
瑞芳鎮、基隆市、汐止鎮、台北市,最後於關渡注入淡水河,與新店溪、大漢溪 同為淡水河三大支流之一。幹流長達89.4 公里,流域面積 490.77 平方公里。從侯 硐介壽橋以下至大華橋為上游段平均坡降約為1/250,以縱谷回春地形著稱,多瀑 布、壺穴、河階景觀,沿線含豐富的礦藏,自七堵大華橋起至南湖大橋為中游段 平均坡降約為 1/4,900,以峽谷、曲流、河階馳名,自南湖大橋到河口為下游段河 床平均坡降約為 1/6,700,河道平緩、故形成帶狀發展的自由曲流地形,穿梭於臺 北盆地之間,為水運及陸運要衝之地。
圖4.1 基隆河流域圖
大致上來說,基隆河主流坡度除瑞芳介壽橋以上之山地河川段較陡外,其餘 坡度均為平緩,支流部分因河川長度較短,且多由山區匯入基隆河,故坡度較陡。
基隆河流域各集水分區面積、河川長度、平均坡度狀況詳見表4.1。基隆河河道分 配圖詳見圖4.2 所示。
表4.1 基隆河流域各集水分區面積、河川長度、平均坡度 主流控制點或
支流分區名稱
流域面積 (平方公里)
河川長度 (公里)
河川平均坡度 備註
介壽橋(瑞芳) 97.83 23.41 1.1/100 山地河川段 大華橋 172.71 43.34 1.7/1000 上游段以上 五堵 180.66 46.53 1.66/1000 中游段以上 南湖大橋 351.15 63.49 1.36/1000 下游段以上 關渡 499.91 83.63 1.05/1000 出海口以上
暖暖溪 25.13 5.17 7/100 支流
大武崙溪 20.74 9.94 5/1000 支流
瑪陵坑溪 18.26 9.36 2.7/100 支流
鹿寮溪 26.10 9.65 4.2/100 支流
保長坑溪 18.66 7.53 6/100 支流
茄苳溪 5.99 3.95 7.7/100 支流
北港溪 17.84 9.24 3.7/100 支流
康誥坑溪 8.98 4.41 6.6/100 支流
叭嗹溪 17.84 8.4 4/100 支流
資料來源:基隆河整體治理計畫(前期計畫)結案報告書
介壽橋
鐵路橋
暖江橋
南湖大橋 百齡橋
關 渡 大 橋
員山子分洪 分 洪
堰 攔河堰
山澗 深澳坑溪 鰈魚坑溪
瑪陵坑溪
北港溪
內溝溪
內雙溪
大坑溪 保長坑溪
暖暖溪
淡 水 河
大華橋
江北橋
社后橋
圖4.2 基隆河河道分配圖
4.2 集水區出水口大華橋之水位預報模式驗證
本研究分析範圍在基隆河上游集水區流域,並依據經濟部水利署第十河川局 公佈之基隆河流域雨量站,利用該流域之四個雨量站,繪出徐昇多邊形網。又基 於出水口設於大華橋水位站,所以只需考量其上游所劃分出徐昇多邊形網之雨量 站,也就是火燒寮、三貂嶺、瑞芳及五堵雨量站。徐昇多邊形網之劃分如圖4.3 所 示。其相關雨量站座標位置詳述如表4.2 所示。
表4.2 雨量站座標位置
雨量站名 UTM-X UTM-Y
五堵 319447.5 2774911
瑞芳 330186.2 2778443
火燒寮 324813.7 2764093
三貂嶺 331962.2 2772822
圖4.3 基隆河大華橋集水區徐昇多邊形網
4.2.1 模式率定
吳等在2004 發表其遞迴函數式:
0
1 1 0
( ) ( ) ( )
N M N
i ij i
i i j
L t a a L t i t b
ξt j t
= = =
= + ∑ − Δ + ∑∑ − Δ
(4-1) 其中,a 、0 a 、i b 為待率定的係數。ijN
為6,M 為 3(三個雨量站:五堵、瑞 芳、火燒寮)。表4.3 為採用公式(4-1),運用最小平方法所求得出的係數值,各個係數詳細資 料如表4.3 所示。
表4.3 率定之係數值
a 6 0.06661 b 12 0.00346 b 36 0.00126 a 5 -0.33171 b 11 0.01022 b 35 -0.00714 a 4 0.17233 b 10 0.00508 b 34 0.00654 a 3 -0.02836 b 26 0.00195 b 33 0.00235 a 2 -0.10951 b 25 0.00496 b 32 -0.00256 a 1 1.17601 b 24 -0.00351 b 31 0.00137 b 16 -0.00256 b 23 -0.00089 b 30 0.00080 b 15 -0.00569 b 22 0.00166 a 0 0.54192 b 14 -0.00816 b 21 0.00751
b 13 0.00165 b 20 0.00171
本文改建其遞迴函數式為:
0
1 1 0 0
( ) ( ) ( ) ( )
N M N N
i ij i k
i i j k
L t a a L t i t b
ξt j t c h t k t
= = = =
= + ∑ − Δ + ∑∑ − Δ + ∑ − Δ
(3-5) 其中,a 、0 a 、i b 、ij c 為待率定的係數。k M 為 4(四個雨量站:五堵、瑞芳、三貂嶺、火燒寮)。
表4.4 為採用公式(3-5),運用最小平方法所求得出的係數值,各個係數詳細資 料如表4.4 所示。
表4.4 率定之係數值
a 6 -0.07159 b 25 -0.00103 b 45 0.00084 a 5 0.13189 b 24 -0.00058 b 44 0.00200 a 4 -0.15274 b 23 -0.00027 b 43 -0.00227 a 3 0.30383 b 22 -0.00523 b 42 -0.00027 a 2 -0.59385 b 21 -0.00200 b 41 0.00538 a 1 1.32254 b 20 0.01670 b 40 -0.00586 b 16 -0.00277 b 36 0.00067 c 6 -0.10270 b 15 -0.00083 b 35 -0.00171 c 5 -0.10603 b 14 0.00000 b 34 0.00011 c 4 0.22087 b 13 -0.00064 b 33 0.00140 c 3 0.02531 b 12 -0.00231 b 32 0.00283 c 2 -0.83563 b 11 0.00439 b 31 -0.00192 c 1 0.96539 b 10 0.00998 b 30 -0.00049 c 0 -0.08620 b 26 -0.00001 b 46 -0.00282 a 0 -4.20301
在率定組部分,圖4.4 為率定組各颱風所預測未來一小時之大華橋水位,藍色 是實測水位,紅色是率定水位,其中舉出辛樂克颱洪事件為例(如圖 4.5 所示),檢 驗率定出的參數計算是否合乎實際觀測水位。圖中虛線為各時刻推算未來一小時 內大華橋水位變化之歷線圖。由結果可知本模式可準確地將降雨量、員山子水位 變化反應到河川水位上,在趨勢及水位高低方面均能掌握。
參數率定
9.00 10.00 11.00 12.00 13.00 14.00 15.00
05/8/4 20:00 05/8/5 20:00 05/9/1 21:00 06/9/10 23:00 07/11/26 10:0007/11/27 10:00 08/9/12 11:00 08/9/13 11:00 08/9/14 11:00 08/9/15 11:00 Time
Stage(M)
實測水位 率定水位
圖4.4 大華橋參數率定
SINLAKU
8.00 9.00 10.00 11.00 12.00 13.00 14.00 15.00
08/9/11 1:00
08/9/11 13:00
08/9/12 1:00
08/9/12 13:00
08/9/13 1:00
08/9/13 13:00
08/9/14 1:00
08/9/14 13:00
08/9/15 1:00
08/9/15 13:00
08/9/16 1:00
08/9/16 13:00
08/9/17 1:00
08/9/17 13:00 Time
Stage(M)
observation with diversion tunnel without diversion tunnel
圖4.5 辛樂克颱洪事件大華橋水位站
表4.5 所示為應用模擬退火演算法,所計算出的參數值。一開始給定之初始值 為使用最小平方法所求出的數值,由表中可發現,經過歷場颱洪事件的發生,其 最佳參數也會有所變動。
表4.5 搜尋最佳參數
優選參數 X( )
α
( )Xβ
( )Xγ
最小平方法 -0.83563 0.96539 -0.08620辛樂克 -0.83560 0.96546 -0.08616
柯羅莎 -0.83572 0.96519 -0.08628
聖帕 -0.83567 0.96530 -0.08623
薔蜜 -0.83554 0.96557 -0.08611
4.2.2 模式驗證
在模式驗證方面,以聖帕、柯羅莎與薔蜜颱洪事件為例,如圖4.6~圖 4.8 所 示。由結果可知利用本模式,一旦由過去歷史記錄找出河川水位與臨前水文條件,
包括降雨量及員山子水位等之間的函數關係後,即可應用於其它的洪水事件,並 準確地將影響水位變化等因素反應到河川水位歷線上。
SEPAT
8.50 9.50 10.50 11.50 12.50
2007/8/17 1:00
2007/8/17 9:00
2007/8/17 17:00
2007/8/18 1:00
2007/8/18 9:00
2007/8/18 17:00
2007/8/19 1:00
2007/8/19 9:00
2007/8/19 17:00
2007/8/20 1:00
2007/8/20 9:00
2007/8/20 17:00 Time
Stage(M)
observation with diversion tunnel without diversion tunnel
圖4.6 聖帕颱洪事件大華橋水位站
KROSA
8.00 9.00 10.00 11.00 12.00 13.00 14.00 15.00 16.00
2007/10/5 1:00
2007/10/5 7:00
2007/10/5 13:00
2007/10/5 19:00
2007/10/6 1:00
2007/10/6 7:00
2007/10/6 13:00
2007/10/6 19:00
2007/10/7 1:00
2007/10/7 7:00
2007/10/7 13:00
2007/10/7 19:00 Time
Stage(M)
observation with diversion tunnel without diversion tunnel
圖4.7 柯羅莎颱洪事件大華橋水位站
JANGMI
8.00 9.00 10.00 11.00 12.00 13.00 14.00 15.00
2008/9/27 14:00
2008/9/27 18:00
2008/9/27 22:00
2008/9/28 02:00
2008/9/28 06:00
2008/9/28 10:00
2008/9/28 14:00
2008/9/28 18:00
2008/9/28 22:00
2008/9/29 2:00
2008/9/29 6:00
2008/9/29 10:00
2008/9/29 14:00 Time
Stage(M)
observation with diversion tunnel without diversion tunnel
圖4.8 薔蜜颱洪事件大華橋水位站
表4.6 大華橋洪水位各颱風預報一小時誤差表
組別 颱洪事件
平均誤差(M) 最大誤差(M)
with diversion
tunnel
without diversion
tunnel
with diversion
tunnel
without diversion
tunnel 率定組 辛樂克 0.065 0.138 0.286 0.642
驗證組
聖帕 0.049 0.074 0.389 0.537
柯羅莎 0.097 0.238 0.279 0.710
薔蜜 0.089 0.185 0.306 0.599
誤差圖
0.00 0.05 0.10 0.15 0.20 0.25
辛樂克 聖帕 柯羅莎 薔蜜
颱風
平均誤差(M )
with diversion tunnel without diversion tunnel
圖4.9 大華橋洪水位各颱風預報一小時誤差
表 4.7 為運用最小二乘法與結合模擬退火演算法互相作比較所計算出的結 果,由圖 4.10 可觀察出最小平方法結合模擬退火演算法,雖與最小平方法計算得 出的結果差異不大,不過卻可較準確地預報大華橋未來一小時之水位。表示可以 將模擬退火演算法加入預報模式中。
表4.7 最小二乘法與結合 SA 比較平均誤差(M)
組別 颱洪事件 最小二乘法(M) 最小二乘法+SA(M)
率定組 辛樂克 0.066 0.065
驗證組
聖帕 0.05 0.049
柯羅莎 0.11 0.097
薔蜜 0.09 0.089
0.066
0.05
0.09
0.065
0.049
0.097
0.089 0.107
0 0.02 0.04 0.06 0.08 0.1 0.12
辛樂克 聖帕 柯羅莎 薔蜜
颱風
平均誤差(M)
最小平方法 最小平方法+SA
圖4.10 最小二乘法與結合 SA 比較平均誤差圖
4.3 淡水河全流域不恆定流河川水位預報模式驗證
淡水河流域位於本島北部,由於都市化的發展,造成其集水區過度的開發利 用,每當颱洪時期,往往有洪患之威脅。淡水河有其三大主要支流,分別為大漢 溪、新店溪與基隆河。三條主要支流中下游地區,為台灣人口最稠密之都會區,
洪水災害一旦發生即有可能造成生命財產上的損失。有鑑於此,建置一套洪水預 報及淹水預警系統,就顯得相形重要。
本文採用的水理分析模式為賴經都(1999)研發的全流域河川不恆定流模式。該 模式是一個實用性強、穩定度高、計算準確的水理模擬模式,可同時處理感潮流、
洪 水 流 及 天 然 或 人 為 導 致 的 急 變 流 , 並 且 可 應 用 於 複 合 - 複 雜 的 渠 道 系 統 (compound-complex channel system)。在起始條件已知的情況下,只要掌握模擬河 段之上下游邊界,該模式即可確實模擬淡水河三大支流(包括感潮河段)洪水運移及 河段中之洪水流況(蔡丁貴與賴經都,2002),其詳細推導如附錄一所示。因此,該 洪水預報工作最主要之重點在準確預測河川洪水位之動態模擬。
該水理分析模式需要其各河川上游邊界點當作輸入條件,現階段上游邊界點 為基隆河的大華橋水位站,新店溪的中正橋水位站,大漢溪的新海橋水位站。由 於基隆河上游員山子分洪道之興建,造成大華橋水位站有高估的現象,故應用本 研究前述之預測大華橋水位站之模式配合水理分析模式,以計算基隆河下游地區 各水位站,以探討本模式之適用性。
考量颱洪時期的中正橋與新海橋水位站需要有其水位資料記錄,以提供水理 分析模式其各河川上游邊界點作為輸入條件,故選擇資料完整且具連續性較高的 聖帕颱風與柯羅莎颱風來進行分析比較。
五堵
4.00 5.00 6.00 7.00 8.00 9.00 10.00 11.00 12.00
2007/8/14 0:00 2007/8/15 0:00 2007/8/16 0:00 2007/8/17 0:00 2007/8/18 0:00 2007/8/19 0:00 2007/8/20 0:00 2007/8/21 0:00 2007/8/22 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.11 聖帕颱風五堵水位站
江北橋
0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.00
2007/8/14 0:00 2007/8/15 0:00 2007/8/16 0:00 2007/8/17 0:00 2007/8/18 0:00 2007/8/19 0:00 2007/8/20 0:00 2007/8/21 0:00 2007/8/22 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.12 聖帕颱風江北橋水位站
社后橋
0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00
2007/8/14 0:00 2007/8/15 0:00 2007/8/16 0:00 2007/8/17 0:00 2007/8/18 0:00 2007/8/19 0:00 2007/8/20 0:00 2007/8/21 0:00 2007/8/22 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.13 聖帕颱風社后橋水位站
表4.8 聖帕颱風於基隆河各水位站預報一小時誤差表
基隆河水位站
with diversion tunnel without diversion tunnel
平均誤差(M) 最大誤差(M) 平均誤差(M) 最大誤差(M)
五堵 0.17 1.15 0.18 1.08
江北橋 0.20 0.9 0.24 0.82
社后橋 0.25 0.96 0.23 1.06
誤差圖
0.00 0.20 0.40 0.60 0.80 1.00
五堵 江北橋 社后橋
水位站
平均誤差(M )
with diversion tunnel without diversion tunnel
圖4.14 聖帕颱風基隆河各水位站洪水位預報一小時誤差
五堵
4.00 5.00 6.00 7.00 8.00 9.00 10.00 11.00 12.00
2007/10/1 0:00
2007/10/2 0:00
2007/10/3 0:00
2007/10/4 0:00
2007/10/5 0:00
2007/10/6 0:00
2007/10/7 0:00
2007/10/8 0:00
2007/10/9 0:00
2007/10/10 0:00
2007/10/11 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.15 柯羅莎颱風五堵水位站比較圖
江北橋
0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00 10.00
2007/10/1 0:00
2007/10/2 0:00
2007/10/3 0:00
2007/10/4 0:00
2007/10/5 0:00
2007/10/6 0:00
2007/10/7 0:00
2007/10/8 0:00
2007/10/9 0:00
2007/10/10 0:00
2007/10/11 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.16 柯羅莎颱風江北橋水位站比較圖
社后橋
0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00 9.00
2007/10/1 0:00
2007/10/2 0:00
2007/10/3 0:00
2007/10/4 0:00
2007/10/5 0:00
2007/10/6 0:00
2007/10/7 0:00
2007/10/8 0:00
2007/10/9 0:00
2007/10/10 0:00
2007/10/11 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.17 柯羅莎颱風社后橋水位站比較圖
南湖大橋
0.00 1.00 2.00 3.00 4.00 5.00 6.00 7.00 8.00
2007/10/1 0:00
2007/10/2 0:00
2007/10/3 0:00
2007/10/4 0:00
2007/10/5 0:00
2007/10/6 0:00
2007/10/7 0:00
2007/10/8 0:00
2007/10/9 0:00
2007/10/10 0:00
2007/10/11 0:00 Time
stage(M)
observation
without diversion tunnel with diversion tunnel
圖4.18 柯羅莎颱風南湖大橋水位站比較圖
表4.9 柯羅莎颱風於基隆河各水位站預報一小時誤差表
基隆河水位站
with diversion tunnel without diversion tunnel
平均誤差(M) 最大誤差(M) 平均誤差(M) 最大誤差(M)
五堵 0.17 0.46 0.27 0.69
江北橋 0.62 1.11 0.83 1.65
社后橋 0.41 0.97 0.58 1.52
南湖大橋 0.54 1.90 0.63 1.45
誤差圖
0.00 0.20 0.40 0.60 0.80 1.00
五堵 江北橋 社后橋 南湖大橋
水位站
平均誤差(M )
with diversion tunnel without diversion tunnel
圖4.19 柯羅莎颱風基隆河各水位站洪水位預報一小時誤差
4.4 員山子分洪流量 4.4.1 員山子分洪工程概述
員山子分洪位於台北縣瑞芳鎮境內,工程起點位於基隆河上游瑞芳鎮瑞柑新 村旁,在主流上築一低型攔河堰,並向東北方闢一條直徑12 公尺長約 2.8 公里之 分洪隧道,出口位置約位於距深澳港東邊約 1.8 公里處匯入東海,下圖為員山子 分洪設施之位置圖。
圖4.20 員山子分洪設施位置(資料來源:經濟部水利署)
員山子分洪之主要工程結構詳細如下(參考資料:經濟部水利署):
(1) 側流堰:為混凝土重力式平頂堰,堰總長 184 公尺,堰高 2.5 公尺,堰頂 處標高 62.5 公尺,採自由溢流方式,溢流後接分洪靜水池。
(2) 攔河堰:屬於混凝土重力式寬頂堰,堰總長 30 公尺,堰高 8 公尺,堰頂 處標高68 公尺。
(3) 分洪靜水池:池底標高為 60 公尺,通往基隆河下游河道設排砂道閘門二 座,通往分洪隧道入口設置分洪堰。
(4) 分洪堰:圓弧型臥箕式堰,堰高 3 公尺,堰頂總長 80 公尺,堰頂處標高 63 公尺,採自由溢流方式,溢流後接束縮段。
(5) 束縮段及入口漸變段:明渠結構,底標高由 60 公尺起以百分之十固定坡 度漸降,接分洪隧道入口。
(6) 分洪隧道:隧道直徑 12m 長約 2.5 公里,分為:
○1 一號分洪隧道:長度1117m,縱坡 1/100。
○2 九份溪過河段:長度27m,縱坡 1/100。
○3 二號分洪隧道:長度1308m,縱坡 1/100。
(8) 出口漸變段及陡槽段:暗渠結構,漸變段底標高以百分之一坡度漸降,陡 槽段底縱坡為一比一,內淨寬由12 公尺漸變至 35 公尺。
(9) 出口消能池:暗渠及明渠結構,內淨寬 35 公尺,內淨高 20.8 公尺;池底 標高負6 公尺、尾檻標高 4 公尺,平時維持 10 公尺水墊以作為消能之用。
圖4.21 員山子分洪設施名稱
(參考資料:基隆河整體治理計畫結案報告)
EL. 63.00 EL. 68.00
EL. 63.00 EL. 68.00
EL. 62.50
N S
E W
LS-1
LS-2 LS-3
LS-4
LS-5
LS-6 LS-7
LS-8
LS-9
LS-X
圖4.22 員山子分洪設施示意圖 (參考資料:亞太儀器有限公司)
表4.10 水壓式水位計
測站 建置點 觀測技術 觀測筒型式 廠牌 儀器型號
LS-1 攔河堰上游左岸 水壓式 鋼筋混凝土管沉箱 KPSI 730-14E-0030 LS-2 攔河堰排洪道 水壓式 PVC 管-預埋方式 KPSI 730-14E-0030 LS-3 攔河堰排砂道 水壓式 PVC 管-預埋方式 KPSI 730-14E-0030 LS-4 分洪堰左側 水壓式 PVC 管-預埋方式 KPSI 730-14E-0030 LS-9 攔河堰下游右岸 水壓式 鋼筋混凝土管沉箱 KPSI 730-14E-0030 (參考資料:亞太儀器有限公司)
表4.11 雷達波式水位計
測站 建置點 觀測技術 觀測筒型式 廠牌 儀器型號
LS-5 進水口隧道洞口 雷達波式 固定架 OTT Kalesto LS-6 進水口束縮段右岸 雷達波式 固定架 E+H FMR 240 LS-7 進水口束縮段右岸 雷達波式 固定架 E+H FMR 240 LS-8 進水口束縮段右岸 雷達波式 固定架 E+H FMR 240 (參考資料:亞太儀器有限公司)