(57)【特許請求の範囲】
【請求項１】
シミュレーション対象領域を基本単位となるメッシュに分割して、前記シミュレーション対象領域をそれぞれのメッシュに割り当てて、シミュレーション対象領域における地形変化のシミュレーションを行う地形変化シミュレーション方法であって、
メッシュ（ｘ，ｙ）毎に対して流下距離に係る流下距離マップＳ（ｘ，ｙ）を定義する流下距離マップ定義ステップと、
前記流下距離マップＳ（ｘ，ｙ）と閾値Ｓslp、Ｓriver、Ｓslp、Ｓseaとに基づいて、メッシュ（ｘ，ｙ）が模擬する地形種別が、斜面域であるか、河川域であるか、斜面域と河川域との間の遷移域であるか、海域であるかを判定する判定ステップと、
メッシュ（ｘ，ｙ）で生産された土砂がどの範囲にまで運搬され堆積するかに係る数値である堆積範囲Ｗを算出する堆積範囲算出ステップと、を有し、
前記判定ステップで、メッシュ（ｘ，ｙ）が河川域であると判定された場合には、Ｗ＝Ｗriver×（Ｓ（ｘ，ｙ）／Ｓriver）αによる計算式で堆積範囲Ｗを計算する（ただし、Wriverは河川先頭での堆積範囲（初期堆積範囲）を示す定数であり、αは粒径係数である）ことを特徴とする地形変化シミュレーション方法。
【請求項２】
前記判定ステップで、メッシュ（ｘ，ｙ）が斜面域であると判定された場合には、Ｗ＝Ｗslpによる計算式で堆積範囲Ｗを計算し（ただし、Ｗslpは定数）、
前記判定ステップで、メッシュ（ｘ，ｙ）が遷移域であると判定された場合には、Ｗ＝｛（１−ｆ）×Ｗslp＋ｆ×Ｗriver｝による計算式で堆積範囲Ｗを計算する（ただし、ｆ＝（Ｓ−Ｓslp）／（Ｓriver−Ｓslp））ことを特徴とする請求項１に記載の地形変化シミュレーション方法。
【請求項３】
メッシュ（ｘ，ｙ）における侵食量マップＶ（ｘ，ｙ）を算出する侵食量マップ算出ステップと、
前記侵食量マップＶ（ｘ，ｙ）を堆積範囲Ｗで除することにより、落水方向における１メッシュあたりの堆積量Ｄｅｐｏを算出する１メッシュあたり堆積量算出ステップと、をさらに有することを特徴とする請求項１又は請求項２に記載の地形変化シミュレーション方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、斜面と河川を含むシミュレーション対象領域における将来の地形変化のシミュレーションを行う地形変化シミュレーション方法に関する。
【背景技術】
【０００２】
地下水などの地質環境の長期的な安定性を評価するためには、地形勾配や土被りが将来的にどの程度まで変化するのかを予測することが重要な課題となる。このような課題を解決するためのアプローチとしては、地形変化のプロセスを明らかにし、一連のプロセスの中で現在がどういう段階にあるかを把握し、さらに将来的に地形がどのように変化していくかを推定することが必要となる。将来の地形変化のうち、河川における侵食、すなわち河床低下は比較的変動が大きいため地下水流動に強く影響するものと考えられており、この変化を推定することは地質環境の安定性評価の重要な１項目となる。
【０００３】
このような観点から、ある程度の広がりを持った地形の大局的な変化を推定する目的のもと、河川と山地斜面を統合した地形変化のシミュレーションが行われてきた。地形変化シミュレーションは、対象領域を基本単位となるメッシュに分割して、これに地形変化のプロセスを模擬する数値計算を行ってメッシュの標高を変化させるものである。具体的には、各メッシュに対し、落水方向の決定、斜面域・河川域などの地形種の決定、侵食量算出、堆積量算出の各ステップからなるシミュレーションを行い、地形変化を模擬する。非特許文献１、２及び３には、河川と山地斜面を統合した地形変化シミュレーションに関連する技術が開示されており、当該非特許文献１及び２では、メッシュの標高値をもとに落水方向、流域面積及び落水勾配を定義する技術（流下距離定義ステップ）及び、地形種別に応じてメッシュ毎の侵食量を算出する技術（侵食土砂量算出ステップ）が記載されている。また、当該非特許文献３では、流域面積に応じてメッシュが模擬する地形種別を判定する技術（判定ステップ）を含め、地質分布を考慮した河床縦断面形（河川の流下距離に対する河川勾配の変化曲線）の地形変化シミュレーション技術が開示されている。非特許文献１、２及び３においては、いずれも侵食による生産土砂量は落水方向において１つ下流側のメッシュに移動し堆積させる（堆積ステップ）ことが記載されている。また、当該非特許文献３では、侵食形態の違いは地質分布で説明可能であり、シミュレーションを用いることで、約百万年前の推定地形面をもとに現在の河川の地形特徴が復元できることが記載されている。
【非特許文献１】野上道男「細密DEMの紹介と流域地形計測」（１９９５），地理学評論，６８（A）−７，ｐ．４６５−４７４
【非特許文献２】野上道男「地理学におけるシミュレーション」（２００５），地理学評論，７８−３，ｐ．１３３−１４６
【非特許文献３】三箇智二、安江健一「河床縦断形のシミュレーション」（２００８），地形，２９（１），ｐ．２７−４９
【発明の概要】
【発明が解決しようとする課題】
【０００４】
非特許文献１、２及び３に記載の技術においては、河川域の地形変化をシミュレーションする際には、一般的に河床縦断面形が指数式で近似し得ることを前提とし、このような指数式近似を成立させるように、土砂移動量の演算式などをたてるようにしていた。
【０００５】
しかしながら、上記のようなシミュレーション方法には、「侵食によって土砂が生産され、生産された土砂が運搬され堆積する」ことを模擬することによって、結果として、河
床縦断面形の形状を得る、という実現象のプロセスをトレースするシミュレーション方法ではなく、河床縦断面形の形状が指数式により近似されるという結果を係数により与えるシミュレーション方法である、という問題があった。また、上記の従来のシミュレーション方法では、河川地形におけるわずかな勾配変化によって大量の土砂移動が出現することなど、シミュレーション結果において実現象とは合致しない結果が得られることがあり、問題となっていた。
【課題を解決するための手段】
【０００６】
上記のような問題点を解決するために、請求項１に係る発明は、シミュレーション対象領域を基本単位となるメッシュに分割して、前記シミュレーション対象領域をそれぞれのメッシュに割り当てて、シミュレーション対象領域における地形変化のシミュレーションを行う地形変化シミュレーション方法であって、メッシュ（ｘ，ｙ）毎に対して流下距離に係る流下距離マップＳ（ｘ，ｙ）を定義する流下距離マップ定義ステップと、前記流下距離マップＳ（ｘ，ｙ）と閾値Ｓslp、Ｓriver、Ｓslp、Ｓseaとに基づいて、メッシュ（ｘ，ｙ）が模擬する地形種別が、斜面域であるか、河川域であるか、斜面域と河川域との間の遷移域であるか、海域であるかを判定する判定ステップと、メッシュ（ｘ，ｙ）で生産された土砂がどの範囲にまで運搬され堆積するかに係る数値である堆積範囲Ｗを算出する堆積範囲算出ステップと、を有し、前記判定ステップで、メッシュ（ｘ，ｙ）が河川域であると判定された場合には、Ｗ＝Ｗriver×（Ｓ（ｘ，ｙ）／Ｓriver）αによる計算式で堆積範囲Ｗを計算する（ただし、Wriverは河川先頭での堆積範囲（初期堆積範囲）を示す定数であり、αは粒径係数である）ことを特徴とする。
また、請求項２に係る発明は、請求項１に記載の地形変化シミュレーション方法において、前記判定ステップで、メッシュ（ｘ，ｙ）が斜面域であると判定された場合には、Ｗ＝Ｗslpによる計算式で堆積範囲Ｗを計算し（ただし、Ｗslpは定数）、前記判定ステップで、メッシュ（ｘ，ｙ）が遷移域であると判定された場合には、Ｗ＝｛（１−ｆ）×Ｗslp＋ｆ×Ｗriver｝による計算式で堆積範囲Ｗを計算する（ただし、ｆ＝（Ｓ−Ｓslp）／（Ｓriver−Ｓslp））ことを特徴とする。
また、請求項３に係る発明は、請求項１又は請求項２に記載の地形変化シミュレーション方法において、メッシュ（ｘ，ｙ）における侵食量マップＶ（ｘ，ｙ）を算出する侵食量マップ算出ステップと、前記侵食量マップＶ（ｘ，ｙ）を堆積範囲Ｗで除することにより、落水方向における１メッシュあたりの堆積量Ｄｅｐｏを算出する１メッシュあたり堆積量算出ステップと、をさらに有することを特徴とする。
【発明の効果】
【０００８】
従来の地形変化シミュレーション方法は、本来シミュレーションの結果として表れるべき河床縦断面形の形状を指数式により近似し指定するというものであったが、本発明の地形変化シミュレーション方法によれば、侵食によって土砂が生産され、生産された土砂が運搬され堆積することを模擬する、という地形変化プロセスにおける実現象と合致したシミュレーションを実施し、その結果として河床縦断面形の形状を得ることが可能となる。また、本発明の地形変化シミュレーション方法によれば、前記の従来のシミュレーション方法において課題となっていた、わずかな勾配による大量の土砂移動などの実現象とは合致しないシミュレーション結果が発生することがない。
【図面の簡単な説明】
【０００９】
【図１】本発明の実施形態に係る地形変化シミュレーション方法を実行させるシステム構成の一例を示す図である。
【図２】本発明の実施形態に係る地形変化シミュレーション方法で用いられるデータベース構成を説明する図である。
【図３】本発明の実施形態に係る地形変化シミュレーション方法における侵食量算出処理のフローチャートを示す図である。
【図４】本発明の実施形態に係る地形変化シミュレーション方法における堆積算出処理のフローチャートを示す図である。
【図５】本発明の実施形態に係る地形変化シミュレーション方法における堆積範囲Ｗの概念を説明する模式図である。
【図６】本発明の実施形態に係る地形変化シミュレーション方法における粒径係数による堆積範囲Ｗの変化を説明する模式図である。
【図７】本発明の実施形態に係る地形変化シミュレーション方法における繰り返し演算を説明する模式図である。
【図８】本発明の実施形態に係る地形変化シミュレーション方法により、非特許文献３の方法によるシミュレーション結果と同様に実際の河床縦断形を再現可能であることを示す図である。
【図９】本発明の実施形態に係る地形変化シミュレーション方法により、氷期における直線的な河床縦断形の形成（実現象）を再現可能であることを示す図である。
【発明を実施するための形態】
【００１０】
以下、本発明の実施の形態を図面を参照しつつ説明する。図１は本発明の実施形態に係る地形変化シミュレーション方法を実行させるシステム構成の一例を示す図である。
【００１１】
図１において、１０はシステムバス、１１はＣＰＵ（Ｃｅｎｔｒａｌ  Ｐｒｏｃｅｓｓｉｎｇ  Ｕｎｉｔ）、１２はＲＡＭ（Ｒａｎｄｏｍ  Ａｃｃｅｓｓ  Ｍｅｍｏｒｙ）、１３はＲＯＭ（Ｒｅａｄ  Ｏｎｌｙ  Ｍｅｍｏｒｙ）、１４は外部情報機器との通信を司る通信制御部、１５はキーボードコントローラなどの入力制御部、１６はディスプレイコントローラなどの出力制御部、１７は外部記憶装置制御部、１８はキーボード、ポインティングデバイス、マウスなどの入力機器からなる入力部、１９はＬＣＤディスプレイなどの表示装置や印刷装置からなる出力部、２０はＨＤＤ（Ｈａｒｄ  Ｄｉｓｋ  Ｄｒｉｖｅ）等の外部記憶装置である。
【００１２】
図１において、ＣＰＵ１１は、ＲＯＭ１３内のプログラム用ＲＯＭ、或いは、大容量の外部記憶装置２０に記憶されたプログラム等に応じて、外部機器と通信することでデータを検索・取得したり、また、図形、イメージ、文字、表等が混在した出力データの処理を実行したり、更に、外部記憶装置２０に格納されているデータベースの管理を実行したり、などといった演算処理を行うものである。
【００１３】
また、ＣＰＵ１１は、システムバス１０に接続される各デバイスを統括的に制御する。ＲＯＭ１３内のプログラム用ＲＯＭあるいは外部記憶装置２０には、ＣＰＵ１１の制御用の基本プログラムであるオペレーティングシステムプログラム（以下ＯＳ）等が記憶されている。また、ＲＯＭ１３あるいは外部記憶装置２０には出力データ処理等を行う際に使用される各種データが記憶されている。ＲＡＭ１２は、ＣＰＵ１１の主メモリ、ワークエリア等として機能する。
【００１４】
入力制御部１５は、キーボードや不図示のポインティングデバイスからの入力部１８を制御する。また、出力制御部１６は、ＬＣＤディスプレイ等の表示装置やプリンタなどの印刷装置の出力制御を行う。
【００１５】
外部記憶装置制御部１７は、ブートプログラム、各種のアプリケーション、フォントデータ、ユーザファイル、編集ファイル、プリンタドライバ等を記憶するＨDＤ（Ｈａｒｄ
Ｄｉｓｋ  Ｄｒｉｖｅ）や、或いは場合によってはフレキシブルディスク（ＦＤ）等の外部記憶装置２０へのアクセスを制御する。
【００１６】
また、通信制御部１４は、ネットワークを介して、外部機器と通信を制御するものであり、これによりシステムが必要とするデータを、インターネットやイントラネット上の外部機器が保有するデータベースから取得したり、外部機器に情報を送信したりすることができるように構成される。
【００１７】
外部記憶装置２０には、ＣＰＵ１１の制御プログラムであるオペレーティングシステムプログラム（以下ＯＳ）以外に、本発明の地形変化シミュレーション方法をＣＰＵ１１上で動作させるシステムプログラム、及びこのシステムプログラムで用いるデータなどがインストールされ保存・記憶されている。なお、上記プログラム用いるデータは外部記憶装置２０に保存しておくこともできるが、場合によっては、このようなデータは通信制御部１４を介してインターネットやイントラネット上の外部機器から取得するように構成することも可能である。
【００１８】
次に、本発明の地形変化シミュレーション方法を実現するシステムプログラムで利用されるデータ構造について説明する。図２は本発明の実施形態に係る地形変化シミュレーション方法で用いられるデータベース構成を説明する図である。本実施形態に係る地形変化
シミュレーション方法では、図２に示すように模擬する地形（シミュレーション対象領域）を基本単位となるメッシュ単位に分割して、シミュレーション対象領域をそれぞれのメッシュ単位に割り当てる。そして、分割された基本単位であるメッシュ単位それぞれを変数とするデータベースを定義する。なお、図２に示す例ではメッシュを矩形形状としたものを示しているが、メッシュ形状は矩形形状に限らず、三角形状などとすることもできる。
【００１９】
本実施形態における地形変化シミュレーション方法では、シミュレーション上、上記のようなメッシュ単位毎に、「河川域」、「斜面域」、これらの中間の「遷移域」、「海域」などの地形種別を割り当て、それぞれのメッシュ単位毎で（模擬上の）「侵食」が発生したり、或いはそれぞれのメッシュ単位毎に（模擬上で）「侵食」によって生産された土砂が（模擬上で）「堆積」したりすることを想定し、このような想定に基づいて地形がどのように変化するかを模擬するものである。なお、本実施形態においては、シミュレーション上における「侵食」、「堆積」などの地形変化模擬現象を、単に「侵食」、「堆積」などの用語で称することがあるが、これらの地形変化現象は言うまでもなく、本実施形態ではあくまでシミュレーション上のものであることを付記しておく。
【００２０】
上記のように定義されたデータベースは外部記憶装置２０に記憶され、地形変化シミュレーション方法を実現するプログラムの処理に基づいて、ＣＰＵ１１によって演算され、適宜書き換えられたりすることによって、各種記憶値の更新が行われる。外部記憶装置２０に記憶されるメッシュ位置（ｘ，ｙ）に依存するデータベース（本実施形態においてはマップなどとも称する）は、図２に示すように、地表標高マップ２０１、落水線ベクトルマップ２０２、流下距離マップ２０３、落水勾配マップ２０４、沖積層厚マップ２０５、岩盤硬度係数マップ２０６、侵食量マップ３００、堆積量マップ４００の８つが定義されている。このうち初期条件として定義されているのは地表標高マップ２０１、沖積層厚マップ２０５、岩盤硬度係数マップ２０６の３つであり、これらの初期条件に基づいて、所定期間にわたる地形の侵食による土砂生産量を計算し、これを侵食マップ３００に記憶・保存し、さらに侵食の過程で発生した土砂の堆積を計算し、これを堆積量マップ４００に記憶・保存する。地形の侵食、土砂の堆積の計算過程で、地表標高マップ２０１、落水線ベクトルマップ２０２、流下距離マップ２０３、落水勾配マップ２０４、沖積層厚マップ２０５の各マップの記憶値が適宜更新される。
【００２１】
上記のマップのうち、地表標高マップ２０１からの落水線ベクトルマップ２０２、流下距離マップ２０３、落水勾配マップ２０４の各マップの作成にかかる考え方については、非特許文献１及び２に基づいている。また、沖積層厚マップ２０５及び岩盤硬度係数マップ２０６の作成にかかる考え方については、非特許文献３の河床縦断形シミュレーションでの「地質係数」の設定をメッシュ位置（ｘ，ｙ）に依存するデータベース（マップ）に拡張したものとする。
【００２２】
以上のように構成される本発明の実施形態に係る地形変化シミュレーション方法を実行させるシステム構成における処理について図３を参照して説明する。図３は本発明の実施形態に係る地形変化シミュレーション方法における侵食量算出処理のフローチャートを示す図である。本フローチャートにて示す侵食量算出処理は、メッシュ（ｘ，ｙ）の１点のみに係るものである。したがって、メッシュで模擬する全てのシミュレーション対象領域について侵食量を算出するには、定義されている全ての（ｘ，ｙ）について同フローチャートによる算出を行うようにする。
【００２３】
ステップＳ１００で、メッシュ（ｘ，ｙ）の侵食量算出処理が開始されると、ステップＳ１０１では、地表標高マップ２０１が参照され、地表標高ｈ（ｘ，ｙ）が読み出される。続いてステップＳ１０２では、この地表標高ｈ（ｘ，ｙ）の値に応じ、標高の低下量が
最大となる隣接メッシュへの方向として、当該メッシュ（ｘ，ｙ）における落水線ベクトルマップ（δｘ，δｙ）（落水線ベクトルマップ２０２）が得られる（非特許文献１に基づく）。
【００２４】
ステップＳ１０３では、落水線ベクトルマップ２０２が参照され、落水線ベクトル（δｘ，δｙ）が読み出される。続いてステップＳ１０４では、この落水線ベクトルの方向に応じ、その流域面積（流入メッシュ数より計算）の平方根より、当該メッシュ（ｘ，ｙ）における流下距離マップＳ（ｘ，ｙ）（流下距離マップ２０３）が得られる（非特許文献３に基づく）。
【００２５】
ステップＳ１０５では、地表標高マップ２０１及び落水線ベクトルマップ２０２が参照され、地表標高ｈ（ｘ，ｙ）及び落水線ベクトル（δｘ，δｙ）が読み出される。続いてステップＳ１０６では、落水線ベクトルに沿う地表標高の変化量δｈより、当該メッシュ（ｘ，ｙ）における落水勾配マップδｈ／δ（ｘ，ｙ）（落水勾配マップ２０４）が得られる（非特許文献１に基づく）。
【００２６】
ステップＳ１０７では、流下距離マップ２０３が参照されて、流下距離Ｓ（ｘ，ｙ）が読み出される。続いてステップＳ１０８では、この流下距離Ｓ（ｘ，ｙ）の値に応じて、当該メッシュ（ｘ，ｙ）が模擬する地形種別が判定される（非特許文献３に基づく）。より詳しくは、このステップＳ１０８では、当該メッシュ（ｘ，ｙ）が模擬する地形種別が、斜面域であるか、河川域であるか、斜面域と河川域との間の遷移域であるか、海域であるかを判定する。
【００２７】
０≦Ｓ（ｘ，ｙ）＜Ｓslpである場合（ステップＳ１０９）には、メッシュ（ｘ，ｙ）
は斜面域であると判定し、Ｓslp≦Ｓ（ｘ，ｙ）＜Ｓriverである場合（ステップＳ１１０）には、メッシュ（ｘ，ｙ）は遷移域であると判定し、Ｓriver≦Ｓ（ｘ，ｙ）である場
合（ステップＳ１１１）には、メッシュ（ｘ，ｙ）は河川域であると判定し、Ｓ（ｘ，ｙ）＝Ｓseaである場合（ステップＳ１１２）には、メッシュ（ｘ，ｙ）は海域であると判
定する。
【００２８】
次の海域以外を模擬するメッシュ（ｘ，ｙ）では、ステップＳ１１４乃至ステップＳ１１６において、メッシュが模擬する地形種別に応じて基準侵食フラックスＰを算出する。ここで、基準侵食フラックスＰについて定義しておく。本実施形態においては、基準侵食フラックスＰが１であるときには、沖積層であれば厚さ１の沖積層を侵食してその分の土砂を生産することができ、岩盤層であればその岩盤硬度係数Ｒ（ｘ，ｙ）に応じて、厚さＰ／Ｒ（ｘ，ｙ）の岩盤層を侵食することができ、その分の土砂を生産することできるものと定義する。なお、沖積層における岩盤硬度係数ＲはＲ（ｘ，ｙ）＝１であるということもできる。
【００２９】
基準侵食フラックスＰの算出においては、ステップＳ１１３において、落水勾配マップ２０４のδｈ／δ（ｘ，ｙ）が参照されて、メッシュ（ｘ，ｙ）の落水方向の隣接メッシュとの間の傾きが考慮される。本実施形態にかかる地形変化シミュレーション方法では、非特許文献２の手法を拡張し、斜面域だけでなく河川域を模擬するメッシュ（ｘ，ｙ）においても、基準侵食フラックスＰの演算式は、指数式を組み込まない単純な拡散の式とする。
【００３０】
ステップＳ１１４では、メッシュ（ｘ，ｙ）が斜面域の場合であって、このときＰ＝ｋslp×ｋriver×δｈ／δ（ｘ，ｙ）による計算式で基準侵食フラックスＰを計算する。ここで、ｋriverは（河川の）メッシュ間隔あたり拡散係数として定義されるものであり、
ｋslpは斜面係数として定義されるものである。ステップＳ１１４においては、河川のメ
ッシュ間隔あたり拡散係数（ｋriver＜１）と斜面係数（ｋslp＜１）と勾配δｈ／δ（ｘ，ｙ）を乗じることによって、基準侵食フラックスＰを計算している。
【００３１】
また、ステップＳ１１６では、メッシュ（ｘ，ｙ）が河川域の場合であって、このときＰ＝ｋriver×δｈ／δ（ｘ，ｙ）による計算式で基準侵食フラックスＰを計算する。
【００３２】
斜面域と河川域との間の領域である遷移域では、非特許文献３の手法に基づき、メッシュ（ｘ，ｙ）におけるＳ（ｘ，ｙ）値が斜面域又は河川域のいずれに近いかに応じて、斜面域又は河川域のいずれの影響が大きいかを考慮するための変数ｆが導入される。ステップＳ１１５に示すように、メッシュ（ｘ，ｙ）が遷移域の場合には、上記のように定義されるｆによって、Ｐ＝｛（１−ｆ）×ｋslp×ｋriver＋ｆ×ｋriver×δｈ／δ（ｘ，ｙ
）｝なる計算式で基準侵食フラックスＰを計算する。
【００３３】
次に、ステップＳ１１７乃至ステップＳ１２１においては、以上のように算出された基準侵食フラックスに応じて、メッシュ（ｘ，ｙ）における生産土砂量が算出される。
【００３４】
ステップＳ１１７では、ステップＳ１１８の判定のために沖積層厚マップ２０５が参照される。ステップＳ１１８では、ステップＳ１１４乃至ステップＳ１１６のいずれかで算出されたＰと、ステップＳ１１７で参照されたＤalv（ｘ，ｙ）とを比較して、Ｐ≦Ｄalv（ｘ，ｙ）の関係が成立するか否かが判定される。
【００３５】
ステップＳ１１８の判定の結果がＹＥＳであるときには基準侵食フラックスＰによって沖積層のみが侵食されることになる。このとき、ステップＳ１２０に進み、沖積層のみからの生産土砂量が算出される
一方、ステップＳ１１６の判定の結果がＮＯであるときには基準侵食フラックスＰによって沖積層に加え岩盤層も侵食されることになる。このとき、ステップＳ１２１に進み、沖積層からの生産土砂量Ｖalvと、岩盤層を削って生産される土砂量Ｖrockとが算出され
る。生産土砂量Ｖrockを算出する上では、ステップＳ１１９において、岩盤硬度係数マップ２０６が参照される。ステップＳ１２１では、メッシュ（ｘ，ｙ）における生産土砂量はＶalvとＶrockとの和として求められる。
【００３６】
上記のようにメッシュ（ｘ，ｙ）における生産土砂量が算出されると、ステップＳ１２２ではこれに基づいて沖積層厚マップ２０５の更新が行われ、ステップＳ１２３では地表標高マップ２０１の更新が行われる。
【００３７】
以上のようにして、ステップＳ１２４でメッシュ（ｘ，ｙ）における侵食量（生産土砂量）に係る侵食量マップＶ（ｘ，ｙ）（侵食量マップ３００）を得ることができる。
【００３８】
以上のようなメッシュ（ｘ，ｙ）の１点に対する侵食量算出処理フローチャートを、シミュレーション対象領域の全域について行うことによって、シミュレーション対象領域全域の侵食量（生産土砂量）を算出することが可能となる。
【００３９】
次に、上記のような侵食量算出処理フローチャートによって算出されたメッシュ（ｘ，ｙ）における生産土砂がどのような形態で堆積するかをシミュレートする堆積算出処理について説明する。図４は本発明の実施形態に係る地形変化シミュレーション方法における堆積算出処理のフローチャートを示す図である。本フローチャートにて示す堆積算出処理は、メッシュ（ｘ，ｙ）の１点で生産された土砂がどの程度の量、どのメッシュへ運搬され堆積するかを計算するものである。したがって、メッシュで模擬する全てのシミュレーション対象領域における堆積の様子を算出するには、定義されている全ての（ｘ，ｙ）について同フローチャートによる算出を行うようにする。
【００４０】
ステップＳ２００で、メッシュ（ｘ，ｙ）の堆積算出処理が開始されると、流下距離マップ２０３が参照されて、流下距離Ｓ（ｘ，ｙ）が読み出される。続いてステップＳ２０２では、この流下距離Ｓ（ｘ，ｙ）の値に応じて、当該メッシュ（ｘ，ｙ）が模擬する地形種別が判定される。より詳しくは、このステップＳ２０２では、当該メッシュ（ｘ，ｙ）が模擬する地形種別が、斜面域であるか、河川域であるか、斜面域と河川域との間の遷移域であるか、海域であるかを判定する。ここで、このような判定のために与えられる閾値については、先に説明したものと同様であるので説明を省略する。
【００４１】
上記のような閾値を用いて、０≦Ｓ（ｘ，ｙ）＜Ｓslpである場合（ステップＳ２０３
）には、メッシュ（ｘ，ｙ）は斜面域であると判定し、Ｓslp≦Ｓ（ｘ，ｙ）＜Ｓriverである場合（ステップＳ２０４）には、メッシュ（ｘ，ｙ）は遷移域であると判定し、Ｓriver≦Ｓ（ｘ，ｙ）である場合（ステップＳ２０５）には、メッシュ（ｘ，ｙ）は河川域
であると判定し、Ｓ（ｘ，ｙ）＝Ｓseaである場合（ステップＳ２０６）には、メッシュ
（ｘ，ｙ）は海域であると判定する。
【００４２】
次に、ステップＳ２０７乃至ステップＳ２０９において、陸域（斜面域、遷移域及び河川域）をを模擬するメッシュ（ｘ，ｙ）で生産された土砂がどの範囲にまで運搬され堆積するかに係る数値（堆積範囲Ｗ）を算出する。
【００４３】
ステップＳ２０７では、メッシュ（ｘ，ｙ）が斜面域の場合であって、このときには、堆積範囲Ｗについては、Ｗ＝Ｗslp（定数）として一定のものと扱うようにする。ただし
、Wslpは斜面域での堆積範囲である。
【００４４】
ここで、堆積範囲Ｗの概念について図５を参照して説明する。図５は本発明の実施形態に係る地形変化シミュレーション方法における堆積範囲Ｗの概念を説明する模式図である。図５では着目メッシュ（ｘ，ｙ）における堆積範囲Ｗが４である場合を例示している。この場合は、図に示すように落水方向である４つ分のメッシュに着目メッシュ（ｘ，ｙ）で生産された土砂の１／４ずつを堆積させるように堆積を模擬する。なお、堆積範囲Ｗの値が４．３などの小数をとるような場合には、落水方向でみて４つのメッシュに着目メッシュ（ｘ，ｙ）で生産された土砂の１／４．３ずつを堆積させ、さらに５つめのメッシュに残りの土砂（生産土砂の０．３／４．３）を堆積させるように堆積を模擬する。
【００４５】
また、ステップＳ２０９では、メッシュ（ｘ，ｙ）が河川域の場合であって、このときには、Ｗ＝Ｗriver×（Ｓ（ｘ，ｙ）／Ｓriver）αによる計算式で堆積範囲Ｗを計算する。ただし、Wriverは河川先頭での堆積範囲（初期堆積範囲）を示す定数であり、αは粒径係数である。図６は本発明の実施形態に係る地形変化シミュレーション方法における粒径係数による堆積範囲Ｗの変化を説明する模式図である。この図６はシミュレーション対象領域の河川の流れ（落水線ベクトルの方向）に沿う断面を模式的に示している。一般的に、河川域における土砂の粒径は上流側で大きく、下流側では小さくなる傾向がある。また、河川流量は流域面積（流下距離Ｓの自乗）の拡大にともない、下流側で増大する傾向がある。このような土砂の粒径及び河川流量の特性によって、河川域における堆積範囲Ｗは上流側で狭く、下流側では広くなるという特徴がある。粒径係数αは、河川域における上記のような土砂の粒径と河川流量に基づく堆積範囲Ｗの特徴を現すために導入される定数であり、およそ４程度の数値である。このような粒径係数αの効果によって、流下距離Ｓ（ｘ，ｙ）が増大すると、（Ｓ（ｘ，ｙ）／Ｓriver）αの項がべき関数的に増大し、同
じく堆積範囲Ｗがべき関数的に増大する。
【００４６】
また、ステップＳ２０８では、メッシュ（ｘ，ｙ）が斜面域と河川域との間の領域である遷移域の場合であって、このときには、メッシュ（ｘ，ｙ）におけるＳ（ｘ，ｙ）値が
斜面域又は河川域のいずれに近いかに応じて、斜面域又は河川域のいずれの影響が大きいかを考慮するための変数ｆが導入される。このような考え方は、侵食量算出処理のフローチャートに関連し説明したものと同様である。堆積算出処理のフローチャートにおいても、遷移域におけるメッシュ（ｘ，ｙ）に対して、変数ｆはｆ＝（Ｓ−Ｓslp）／（Ｓriver−Ｓslp）として定義し、上記のように定義されるｆによって、Ｗ＝｛（１−ｆ）×Ｗslp＋ｆ×Ｗriver｝なる計算式で堆積範囲Ｗを計算する。
【００４７】
以上のようにステップＳ２０７乃至ステップＳ２０９で堆積範囲Ｗが算出されると、続いてステップＳ２１０で侵食量マップＶ（ｘ，ｙ）が参照されて、メッシュ（ｘ，ｙ）における生産土砂量が取得され、さらにステップＳ２１１で侵食量マップＶ（ｘ，ｙ）が堆積範囲Ｗによって除されることで、メッシュ（ｘ，ｙ）からの落水方向における１メッシュあたりの堆積量Ｄｅｐｏが求められる。
【００４８】
ステップＳ２１２からステップＳ２２０は、上記のようにして求められた１メッシュあたりの堆積量Ｄｅｐｏがどのメッシュに堆積するかを求める処理となる。
【００４９】
ステップＳ２１２では、ｎ＝１がセットされて繰り返し演算のルーチンが開始される。この繰り返し演算のルーチンは、例えば、Ｗ＝４であるときには、ｎ＝１からｎ＝４までの４回繰り返されることで、Ｄｅｐｏ（＝Ｖ（ｘ，ｙ）／４）がどのメッシュに堆積されるかを求めるルーチンである。
【００５０】
なお、例えば、Ｗ＝４．３などの自然数でないときには、繰り返し演算はｎ＝１からｎ＝５までの５回繰り返す。ｎ＝１からｎ＝４までの繰り返し演算で、４つのメッシュにＤｅｐｏ（＝Ｖ（ｘ，ｙ）／４．３）ずつを堆積させ、ｎ＝５の演算で、５つめのメッシュに残りの土砂（生産土砂の（０．３／４．３）×Ｖ（ｘ，ｙ）を堆積させる。
【００５１】
以上のような繰り返し演算のルーチンについてより詳しく説明する。まず、ステップＳ２１３では落水線ベクトルマップ２０２が参照されることで落水方向が求められ、これに応じてＤｅｐｏを堆積させる堆積先のメッシュが選定される。
【００５２】
ステップＳ２１４では、次の堆積先であるＳ（ｘ＋δｘ，ｙ＋δｙ）のメッシュの地形種別の判定が行われる。ステップＳ２１４における判定の結果、堆積先が陸域であるときにはステップＳ２１５に進み、堆積先が海域であるときにはステップＳ２１８に進み、堆積先がシミュレーション対象領域外であるときにはステップＳ２２０に進む。
【００５３】
堆積先が海域であるときに進むステップＳ２１８では、（繰り返し演算の途中であっても）河口に相当するメッシュに残りの土砂を堆積させるように模擬し、ステップＳ２１９で河口堆積量マップＤsea（ｘ，ｙ）を作成する。この河口に堆積された土砂については
、これを海底に堆積させるための別のルーチンによって処理がなされるが、本発明は海面より高い地形部分に焦点を当てるものであるので、本実施形態においてはこの海底への拡散処理ルーチンについての詳細な説明は割愛する。
【００５４】
堆積先がシミュレーション対象領域外であるときに進むステップＳ２２０では、（繰り返し演算の途中であっても）領域外への流出口に相当するメッシュに堆積させた分を除く残りの土砂はなくなったものとしてカウントする。
【００５５】
堆積先が陸域であるときに進むはステップＳ２１５では、既に堆積している土砂に、メッシュ（ｘ，ｙ）で生産された土砂分であるＤｅｐｏを堆積する。ステップＳ２１６では、ｎを１インクリメントして、ステップＳ２１７ではｎ＜Ｗ＋１であるか否かが判定され、繰り返し演算の要否が判定される。
【００５６】
ここで、堆積範囲がＷ＝４であり、かつ、堆積先が全て陸域である場合で、ｎ＝１からｎ＝４までの４回の繰り返し演算が行われることを例にとり、繰り返し演算のルーチンで行われる処理について図７を参照して説明する。図７は本発明の実施形態に係る地形変化シミュレーション方法における繰り返し演算を説明する模式図であり、図７（Ａ）から図７（Ｄ）の順で繰り返し処理が進んでいる。
【００５７】
図７（Ａ）に示すように、メッシュ（ｘ，ｙ）で生産された土砂の１／４分の第１の堆積先が、落水線ベクトルマップ２０２によって求められ、求められた堆積先に土砂の１／４分を堆積させる。次の図７（Ｂ）では、第１の堆積先からさらに落水線ベクトルマップ２０２によって第２の堆積先が求められ、この第２の堆積先に、土砂の１／４分を堆積させる。以降、繰り返し演算ルーチンでは、図７（Ｄ）までこの手順を繰り返す演算処理を行う。
【００５８】
ステップＳ２２１では、メッシュ（ｘ，ｙ）で生産された土砂の堆積状況を示す堆積量マップＤ（ｘ，ｙ）（堆積量マップ４００）が求められる。また、このような堆積量マップＤ（ｘ，ｙ）の変化に伴い、ステップＳ２２２では、沖積層深度マップ２０５の更新が行われ、ステップＳ２２３では地表標高マップ２０１の更新が行われる。
【００５９】
以上の堆積算出処理フローチャートは、メッシュ（ｘ，ｙ）の１点で生産された土砂の堆積状況を算出するものであり、この堆積算出処理フローチャートを、シミュレーション対象領域の全域について行うことによって、シミュレーション対象領域全域の堆積の状況を算出することが可能となる。
【００６０】
この過程を繰り返すことによって、指数式を用いた従来の地形変化シミュレーション結果（非特許文献３）と同様に、約百万年前の推定地形面をもとに現在の河川の地形特徴が再現できる（図８）。また、粒径係数を変化させたシミュレーションにより、氷期における河川の運搬力（粒径の減少傾向）の低下による直線的な河床縦断形の形成を再現することができる（図９）。
【００６１】
以上、本発明の地形変化シミュレーション方法によれば、侵食によって土砂が生産され、生産された土砂が運搬され堆積することを模擬する、という地形変化プロセスにおける実現象と合致したシミュレーションが可能となる。また、本発明の地形変化シミュレーション方法によれば、シミュレーション結果（河床縦断形の形成など）において実現象と合致しない状況が発生することがない。
【符号の説明】
【００６２】
１０・・・システムバス、１１・・・ＣＰＵ（Ｃｅｎｔｒａｌ  Ｐｒｏｃｅｓｓｉｎｇ  Ｕｎｉｔ）、１２・・・ＲＡＭ（Ｒａｎｄｏｍ  Ａｃｃｅｓｓ  Ｍｅｍｏｒｙ）、１３・・・ＲＯＭ（Ｒｅａｄ  Ｏｎｌｙ  Ｍｅｍｏｒｙ）、１４・・・通信制御部、１５・・・入力制御部、１６・・・出力制御部、１７・・・外部記憶装置制御部、１８・・・入力部、１９・・・出力部、２０・・・外部記憶装置、２０１・・・地表標高マップ、２０２・・・落水線ベクトルマップ、２０３・・・流下距離マップ、２０４・・・落水勾配マップ、２０５・・・沖積層厚マップ、２０６・・・岩盤硬度係数マップ、３００・・・侵食量マップ、４００・・・堆積量マップ
【図１】
0005422833.tif 000002
【図２】
0005422833.tif 000003
【図３】
0005422833.tif 000004
【図４】
0005422833.tif 000005
【図５】
0005422833.tif 000006
【図６】
0005422833.tif 000007
【図７】
0005422833.tif 000008
【図８】
0005422833.tif 000009
【図９】
0005422833.tif 000010