河川技術論文集,第18巻,2012年6月

論文

実測水位に基づく分布型流出・平面２次元洪水追跡モデルのパラメータ最適化法
PARAMETER OPTIMIZATION METHOD BASED ON OBSERVED WARTER
SURFACE OF FLOOD FLOWS FOR DISTRIBUTED RAINFALL-RUNOFF AND 2D
FLOOD FLOWS SIMULATION MODEL
重枝未玲1・秋山壽一郎2・野村心平3
Mirei SHIGE-EDA, Juichiro AKIYAMA and Shinpei NOMURA
1正会員

博士(工)

九州工業大学大学院准教授 工学研究院建設社会工学研究系（〒804-8550 北九州市戸畑区仙水町1-1）

2フェロー会員

Ph.D.

九州工業大学大学院教授

3学生会員

九州工業大学大学院工学研究院建設社会工学研究系（同上）工学府建設社会工学専攻博士前期課程（同上）The parameter optimization method based on observed water surface of flood flows for distributedrainfall-runoff and 2D flood flows simulation model was examined. The design of experiments andresponse surface method were used for the optimization method. The rain-fall runoffs were simulated
by a distributed hydrological model using a high-resolution scheme. The flood flows in the main river anda tributary were simulated simultaneously by the 2D unsteady free-surface flow model with treatment ofparallel computation and sink term semi-implicitly. The flooding event in the Onga river basin in July
2003, 2009 and 2010 were simulated by the model. The simulated results were compared with theobserved water level. It shows that the optimization method is a useful tool for parameter estimation..
Key Words : distributed rainfall-runoff model, 2D unsteady free-surface flows simulation model,
parameter optimization, design of experiments, Response surface method, Onga river

１． はじめに

近年，観測史上最大やそれに近いような豪雨による水害が頻発し，甚大な被害が生じている．地球温暖化に伴う気候変動により年最大日雨量が増大する予測1),2)もなされており，今後，観測史上最大の豪雨が発生する可能性はますます高くなると考えられる．
このような自然環境の変化とともに，社会的背景から「できるだけダムにたよらない治水対策」への転換3)も打ち出されており，治水対策は，流域全体で流出の抑制対策を行い河道への雨水流出を遅らせ洪水ピーク流量を低減させることで治水安全度を確保する「流域と一体となった治水対策」へと移行している．超過洪水対策を検討する上で，流域の治水施設が超過洪水に及ぼす影響を知ることは重要であり，そのためには，降雨外力から治水施設によって変化する流出プロセスやその結果によって生じる洪水プロセスを予測できるモデルが不可欠である．
流出プロセスの予測については，土地利用に応じてモデルパラメータを設定する地物指向の分布型流出モデル
が，洪水プロセスの予測については，河道の縦横断面形状，平面形状，樹木繁茂状況や河床材料・状態に応じた樹木群透過係数や粗度係数を設定する平面2次元モデル等が開発6),7)されており，流量ハイドログラフや水位ハイドログラフを適切に予測できることを示している．著
者らは，このようなモデルを用いれば，任意の降雨から水位ハイドログラフの予測が可能となり，モデルの検証やモデルパラメータの設定が流量に比べ観測精度の高い水位に基づき行うことが可能になると考え，分布型流出解析モデルと高度な平面2次元解析モデルを用いた洪水追跡法による「分布型流出・平面2次元洪水追跡モデ
ル」を開発し，同モデルが遠賀川流域を対象に本川と複数の支川を対象に，流域と河道特性を踏まえた上で降雨外力から水位ハイドログラフ等の洪水流の挙動を予測できることを示した8),9),10)．しかし，モデルパラメータの最適化には課題が残った．
本研究は，以上のような背景を踏まえ，流域と河道特性を踏まえた上で，降雨外力から本川および支川の水位ハイドログラフなどの洪水流の予測を行なうことが可能4),5)

な「分布型流出解析・平面2次元洪水追跡モデル」の開発を最終的な目的としたものである．ここでは，同モデルの特徴を踏まえ，従来の流出モデルのように流量ハイドログラフではなく，観測所が多くかつ精度が高い水位ハイドログラフに基づき，分布型流出・平面2次元洪水追跡モデルのモデルパラメータの最適化について検討した．

２．モデルと最適化手法の概要

U E F
 S1  S 2  S 3  0tx y(1)図-1 最適化の手順

地盤高と落水線

森林
荒地
田
農地
都市

土地利用

水位観測所地表面部

地盤高4.5ｍ

河道部
100
降雨量(mm)(1) 分布型流出・平面2次元洪水追跡モデルの概要10)
分布型流出・平面2次元洪水追跡モデル10)は，降雨を入力条件として河道での水位ハイドログラフ等を予測できるモデルである．同モデルは，①セル分布型流出解析モデルと②平面2次元自由表面流モデルで構成され，①では標高や土地利用などの流域特性を考慮した上で流域の雨水の挙動の予測を，②では河道平面形状などの河道特性を考慮した上で洪水追跡を行う．両モデルは，①より求まる河道への流出流量を②の境界あるいは内部境界条件として与えることで接続される．
セル分布型流出解析モデルは，斜面流と河道流の洪水追跡をkinematic wave 法で行うモデルである．なお，河道流の追跡は，平面2次元洪水追跡の対象としない河道についてのみ実施する．セル分布型の基礎方程式を式(1)に示す．
h q(1)rt x
ここに，t=時間，x=斜面の距離，h=水深，r=降雨強度，
q=単位幅流量である．単位幅流量qには飽和流，不飽和流，表面流を区別できる斜面表層の流量流積関係式11)を用いた．同モデルのモデルパラメータは，N=等価粗度係数，ka=重力水が卓越するA層内の透水係数，dc＝マトリックス部の最大保水量を水深で表した値，ds＝重力水を含めて表層土中に保水しうる最大水深，=重力水部と不飽和水部との飽和透水係数の比である．式(1)の離散化は，高解像度風上解法の一つであるMUSCL法12)に基づき行われている．
平面2次元洪水追跡モデルは，複数の河川での洪水追跡をDynamic wave法で行うモデルである．同モデルは，有限体積法，流束差分離法(FDS法)，領域分割法による並列計算などのCFD技術によるPSA-FUF-2DFモデル(AParallelized Spatial Averaged Finite volume method onUnstructured grid using FDS technique for 2D Free-surfaceflows)10)である．なお，同モデルの摩擦勾配の取り扱いについては半陰解法を用いている．PSA-FUF-2DFモデルの基礎方程式は，以下に示す通りである．

5000

612時刻(h)1824

図-2 仮想流域と降雨ハイエトグラフ



TU  h, uh, vh  ； E  uh, u 2 h  1 / 2 gh 2 , uvh ;

T



T

F  vh, uvh, v 2 h  1 / 2 gh 2 ；S1  0, ghS ox , ghS oy  ； S 2  0, ghS fx , ghS fy  ；TTS 3  qr , Fx , Fy T

ここに，U=保存量ベクトル，E，F=x，y方向の流束ベク
トル，S1=河床勾配ベクトル，S2=摩擦勾配ベクトル，
S3=発生項・消滅項ベクトル，h = 水深，u，v = x，y方向の流速，g＝重力加速度，qr=単位面積当りの流入流量(－：流入，+：流出)，Sox，Soy = x，y方向の河床勾配，Sfx，Sfy = x，y方向の摩擦勾配，Fx，Fy = 計算メッシュ内に樹木などの物体群が含まれる場合に付加されるx，y方向の流体力項である．摩擦勾配は，Manning の公式を用いて，流体力項は抵抗係数Cdを用いた式13)で計算される．
式(1)の離散化は，時間積分についてはEulerの陽解法に基づき，空間積分については有限体積法と流束差分離法12)に基づき行われている．
分布型流出・平面2次元洪水追跡モデルのモデルパラ

表-1 モデルパラメータの設定範囲

表-2 観測結果として用いたモデルパラメータ-1/3
d s (m) d c (m)k a (m/s)
βN (m ・s)
森林 0.4～0.6 0～0.3 0～0.3 3.5×10-5～0.03 4～6
荒地 0.2～0.4
田
1～30.00
農地 0.2～0.4 0.00
都市 0.01～0.1
水域
0.035

N (m-1/3・s) d s (m) d c (m) k a (m/s) β
森林 0.600.20
0.15
0.015 5.00
荒地 0.30田
2.00.00
0.00
農地 0.30
都市 0.055
水域 0.035

表-3 仮想流域での解析条件の一例-1/3

Run1Run2

森林
0.4
0.8

N (m ・s)
荒地
田
0.31
30.24
1

農地
0.2
0.27

都市
0.1
0.0145d s (m)d c (m)00k a (m/s)
0 0.028502
0 0.000035

β5.8
5.4・・・・・・・・・・・・・・・・・・・・・・・・・・・・・・

Run37
Run38

0.8
0.4

0.33
0.2

1.8
30.2
0.26

0.01
0.1

0.15
0.15

0.15 0.024007
0.3 0.009025

4.7
4・・・

・・・

・・・

・・・

・・・

・・・

・・・

・・・・・・・・・

Run71
Run72

0.4
0.8

0.2
0.4

31

0.4
0.2

0.01
0.1

000.3
0.03
0.3 0.000035

46

イドログラフを予測できるためである．
メータは，流域の土地利用(森林，荒地，田，農地，都最適化の概要は図-1に示す通りであり，その手順は市)に応じた等価粗度Nと，飽和透水係数の比，土層中の大空隙部の飽和透水係数ka，マトリックス部の層厚dc， (1) 実験計画法による解析条件の設定，(2) 解析条件による分布型流出・平面2次元洪水追跡の実行と誤差のカ全空隙に相当する層厚dsと河道の粗度係数nである．な
イ2乗基準の算定，(3) 解析結果に基づく応答曲面の作成
お，，ka，dc，dsは土地利用が森林の場合についてのみと目的関数を最小とするパラメータの推定である．
考慮した．
(2) 最適化手法の概要
分布型流出・平面2次元洪水追跡モデルは，前述した
ように，降雨外力から水位ハイドログラフ等の洪水流の挙動を予測できることが特徴である．ここでは，この特徴を踏まえ，従来の流出モデルのような流量ハイドログラフではなく，観測所が多くかつ精度が高い水位ハイドログラフに基づくモデルパラメータの最適化手法について述べる．
最適化手法には，LM(Levenberg-Marqardt)法14)や群知能最適化手法15)等があり，分布型流出解析モデルの最適なパラメータ値の推定に成功している．本研究では，車両の設計分野の寸法最適化16)で用いられている実験計画法と応答曲面法を適用した．なお，著者らの知る限り，同様な手法の流出モデルへの適用例は存在しない．
実験計画法は，モデルパラメータが目的関数に及ぼす影響を把握するための効率的な解析条件を設定する方法である．応答曲面法は，モデルパラメータと目的関数との関係を応答曲面と呼ばれる2次の多項式で求め，その関係から最適化パラメータを推定する方法である．ここでは，目的関数を観測水位と解析水位を用いたカイ2乗基準とし，最適化を行うモデルパラメータは，5つの土地利用別のN，，ka，dc，dsの9つとした．なお，河道の粗度係数は最適化を行うパラメータには含めていない．
これは，河道の縦横断面形状，平面形状および非定常性に起因して起こる貯留現象を取り扱うことができる平面
2次元モデルであれば，樹木繁茂状況や河床材料や状態に応じた適切なパラメータ値を設定することで，水位ハ３．最適化手法の妥当性の検討
はじめに，図-2に示す仮想的な流域と河道を対象に，仮定したパラメータ値より得られる解析結果を観測水位として，本手法の最適パラメータ値の妥当について検討した．
解析対象流域の流域面積は5.2km2で流域の地形と土地利用は図-2に示す通りである．河道は矩形断面で延長2300ｍ，川幅25ｍ，勾配を9/4600，粗度係数n=0.035とした．降雨外力には図-2に示す後方集中型の降雨を与えた．
モデルパラメータの範囲を表-1に示すように設定し，表-2に示すモデルパラメータによる解析結果を観測値として用いた．なお，観測地点は図-2中の○で示す河道中央とした．
表-3は，表-1のパラメータ値の範囲から，実験計画法に基づき設定した解析条件の一例を示したものである．
この条件で分布型流出解析・平面2次元洪水追跡を行い，応答曲面を算定した．
表-4は，応答曲面より得られる最適パラメータ値を示したものである．これより，(1) マトリックス部の層厚dc，全空隙に相当する層厚dsで15～20％程度の誤差があるが，(2) 本最適化手法は比較的良い精度でモデルパラメータを再現できることがわかる．図-3は，最適パラメータにより求めた水位ハイドログラフの解析結果を示したものである．図中には他のパラメータを用いた結果も示している．これより，(1) 最適パラメータ値を使用

表-4 仮想流域でのモデルパラメータの最適値-1/3
N (m ・s) d s (m) d c (m) k a (m/s) β
森林 0.620.23
0.120.015 5.11
荒地 0.30
田
1.87
0.00
0.00
農地 0.28
都市 0.053
水域 0.035
4最適解
観測値(仮定)
43.5

解析結果

水位(m)

3.5

3

32.5

22.5

206012時刻(h)6121824

図-3 水位ハイドログラフの比較

図-4 遠賀川流域の概要と解析対象区間

図-5 降雨量の空間分布

していない場合，水位上昇時やピーク値や減水時での誤差が生じ，特に減水時に大きいこと，(2) 最適パラメー
タ値を用いた場合，減水時で若干の誤差があるものの，増水時，ピーク時の水位ともに再現性高いことが確認できる．
このように，最適パラメータ値を用いることで十分な精度で水位ハイドログラフを予測可能であることから，本最適化手法は妥当であることがわかる．ただし，減水時に誤差があることから土層水の流れに関するパラメー
タには検討の余知が残る．この点については，今後検討したいと考えている．

４．遠賀川流域への適用
次に，遠賀川流域を対象に観測水位に基づき最適化パラメータを推定するとともに，複数の豪雨に対する同モデルの予測精度と最適化パラメータの妥当性について検討した．

(1) 遠賀川と流域の概要
遠賀川は幹川流路延長61km，流域面積1,026km2の一級河川であり，流域の土地利用は，森林が約52%，水田や果樹園等の農地が約19%，宅地等市街地が約18%，荒地が約6%を占める．図-4に遠賀川流域の航空写真と雨量・水位観測所を示す．
(2) 解析条件の概要
分布型流出解析では遠賀川流域全体を，平面2次元洪水追跡では図-4に示す勘六橋観測所より上流を，それぞれ解析対象とした．
降雨には2003年(パターン①)，2009年(パターン②)，
2010年(パターン③)の3つの降雨パターンを用い，パターン②の降雨でパラメータの最適化を行った後，そのパラメータ値を用いてパターン①と③での解析を行った．
図-5は解析に用いた総降雨量分布を示したものである．
いずれも中流域を中心とした降雨分布であることが特徴である．
分布型流出解析の開始時刻は，それぞれパターン①である．これより， (1) いずれの観測所についても，若干のずれがあるものの，水位ハイドログラフの波形やピーク値などを再現できること，(2) 特に，パターン②と③では2山の波形を持つ複雑な水位ハイドログラフとなるが，これらも再現していること，などがわかる．このよ
うに，本最適化手法は，モデルパラメータを十分な精度は2003年7月18日2時，パターン②では2009年7月24日1時， で推定できることがわかる．しかしながら，(3) パター
ン①では減水が早いこと，(4) パターン②，③では2波目
パターン③では2010年7月12日1時とした．平面2次元洪のピークが早く減水も早いことなどがわかる．これは，
水追跡では，図-4に示す勘六橋観測所より上流の遠賀川
先述の仮想領域での検討結果を踏まえると，本最適化手本川と支川の穂波川を対象とし，河道を12区間に分割し法では土壌水のパラメータの推定が十分ではなく土層厚た．河道横断面図はパターン①では2003年のものを，パなどのパラメータが現実とは異なる可能性があり，表面ターン②と③では2009，2010年のものを用いた．総メッ流が支配的となる増水時に比べ，相対的に土壌水の影響シュ数はパターン①で7743個，パターン②と③で5893個であった．平面2次元洪水追跡のパラメータについては， が大きくなる減水時に差が生じたためと考えられる．
以上から，本最適化手法は，土壌水の流れに関するパ粗度係数には計画で用いられた値0.035(m-1/3s)を，下流端ラメータの推定には課題が残るが，最適パラメータ値をには勘六橋観測所の実測水位を与えた．平面2次元洪水推定する上で有用であることが確認された．
追跡の解析開始時刻は，それぞれパターン①では2003年7月18日15時，パターン②では2009年7月24日12時，パターン③では2010年7月13日1時とした．なお，流出解析
５．おわりにおよび平面2次元洪水追跡の解析データの作成方法の詳
細は，参考文献8),9),10)に示す通りである．実験計画法に基
本研究では，水位ハイドログラフに基づき分布型流づき設定した解析条件は，表-3と同様である．
出・平面2次元洪水追跡モデルのモデルパラメータの最適化について検討した．その結果，以下のことが確認さ(3) 結果と考察れた．
表-5は，分布型流出・平面2次元洪水追跡の解析結果(1) 実験計画法と応答曲面法を用いた実測水位に基づくに基づき求めた応答曲面より得られる最適パラメータ値
最適化手法は，モデルパラメータを十分な精度で推を示したものである．これより，(1) 等価粗度係数につ定できること．
いては，市街地で若干大きな値を取るが，いずれも標準(2) 分布型流出・平面2次元洪水追跡モデルは，最適化値程度であること，(2) マトリックス部の層厚dcと全空パラメータを用いることで，複数の豪雨に対して本隙に相当する層厚dsには大きな差がないこと，などわか川および支川の水位ハイドログラフを十分な精度でる．パターン②の降雨では，小河原雨量観測所等で時間再現できること．
降雨100mmを越えるなどの豪雨であったことを踏まえる(3) ただし，土壌水の流れに関するパラメータの推定にと，森林の保水能力を超えると早期に地表面流が生じ，は課題が残るため，実験計画法による解析条件の追流域からの流出は地表面流が支配的であったと考えられ加など，その改善方法について検討が必要なこと．
る．
今後は，上記のパラメータ推定の改善を行うとともに，図-6は，遠賀川・穂波川について，パターン②の水面地物や土壌の湿潤状況と最適パラメータとを関連付ける
形の経時変化の解析結果と実測値との比較を行ったもの
ことで，土地利用や地質等からパラメータを設定できるである．これらより，本川の遠賀川では，(1) 24日17:00
手法を検討したいと考えている．
～23:00では，勘六橋から穂波川との合流部付近までの区間で水位上昇が顕著となること，(2)一方で，24日
謝辞：本研究を実施するに当たり，国土交通省遠賀川河
23:00～25日20:00では，穂波川合流部よりも上流で水位
川事務所からデータ提供などで，当時本学学部4年生のが上昇した後，低下すること，(3) 解析結果は実測値を網中達郎君，大谷昌成君にはデータ整理等で協力を得た．
概ね再現していること，などが確認できる．また，支川
ここに記して感謝の意を表します．の穂波川では(1) 本川の水位上昇により，合流部から秋松橋水位観測所の区間で水位が上昇すること，(2) 水面
参考文献
形は，増水時，減水時ともに概ね相似形を保ち変化する
1) 国土交通省：地球温暖化に伴う気候変動について，大規模こと，(3) 解析結果は実測値を概ね再現していること，水害対策に関する専門調査会，2007.などが確認できる．
2) 小林健一郎，寶 馨，中北英一：全球気候モデル出力を用
図-7は，パターン①～③について，最適パラメータ値
いた日本域の100年確率日降水量の将来予測，水工学論文を用いた水位ハイドログラフの解析結果を示したもので表-5 遠賀川流域でのモデルパラメータの最適値-1/3
N (m ・s) d s (m) d c (m) k a (m/s) β
森林 0.63 0.133 0.127
0.015 5.09
荒地 0.30
田2.53
0.00
0.00
農地 0.29
都市 0.055
水域 0.035

図-6 水面形の経時変化(パターン② 左：24日17:00～23:00，右：24日23:00～25日20:00)

図-7 各観測所での水位ハイドログラフの比較(左：パターン①，中：パターン②，右：パターン③)集，第54巻，pp.223-228，2010．
3) 今後の治水対策のあり方に関する有識者会議：今後の治水対策のあり方について(中間とりまとめ)，国土交通省河川局，2010．
4) 天口英雄，河村明，高崎忠勝：地物データGISを用いた新たな地物指向分布型都市洪水流出解析モデルの提案，土木学会論文集，Vol.63，No.3，pp.206-223，2007．

追跡と改修効果の評価，土木学会論文集B1(水工学)，
Vol.68，No.4，I_1429-I_1434 ，2012．
11) 立川康人，永谷 言，寶 馨：飽和・不飽和流れの機構を導入した流量流積関係式の開発，水工学論文集，第48巻，pp.7-12，2004．
12) 藤井孝考藏：流体力学の数値計算法，東京大学出版会，
1994．

5) 天口英雄，河村明，荒木千博，高崎忠勝：神田川上流域の

13) 重枝未玲，朝位孝二，坂本 洋，長太茂樹，秋山壽一郎，

地物データを用いた洪水流出モデルの構築とその適用，河樋口直樹，重岡広美，徳永智宏：大野川とその派川の乙津川技術論文集，第15巻，pp.377-382，2009．
6) 福岡捷二，渡邊明英，原俊彦，秋山正人：水面 形の時間変化と非定常二次元解析を用いた洪水流量ハイドログラフと貯留量の高精度推算，土木学会論文集，No.761//II-67，
pp.45-56，2004．
7) 福岡捷二，佐藤宏明，藤澤寛，大沼史佳：洪水流と河道の川を包括した平面2次元洪水流解析と河道内樹木が乙津川への分流量に及ぼす影響，河川技術論文集，第12巻，pp.85-90，2006．
14) 小林健一郎，寶 馨，立川康人：最適化手法による分布型
降雨流出モデルのパラメータ推定，水工学論文集，第51巻，pp.409-414，2007．

樹木繁茂形態に基づく樹木群透過係数と粗度係数の算定法， 15) 小槻峻司，田中賢治，小尻利治，浜口俊雄： 群知能最適
水工学論文集，第51巻，pp.607-612，2007．
8) 重枝未玲，秋山壽一郎，小園裕司，上川謙：遠賀川流域の分布型流出解析と平面2次元洪水追跡，水工学論文集，第54巻，pp.517-522，2010．
9) 重枝未玲，秋山壽一郎，小園裕司：遠賀川流域の分布型流出解析と複数の支川を考慮した平面2次元洪水追跡，河川化手法を用いた分布型流出モデルのパラメーター同定，土木学会論文集B1(水工学)，Vol.68，No.4，I_523-I_528，
2012．
16) 尼子龍幸：ファーストオーダーアナリシス FOAにおける応答曲面法による寸法最適化，豊田中央研究所R&Dレビュー，第37巻，1号，pp.31-36，2002．

技術論文集，第16巻，pp.443-448，2010．
10) 重枝未玲，秋山壽一郎，草野浩之，野村心平：高解像度風上解法を用いた遠賀川流域の分布型流出・平面2次元洪水（2012.4.5受付）