多領域の接続計算法，（2）外洋域の無反射透過スキーム の構築，（3）実地形データに基づく計算の実施とその結 果の検証，の3点について検討を実施する．
2. 格子ボルツマン法
（1）格子ボルツマン法の概要
LBMは分子動力学における格子気体法（Lattice Gas
Automata）から発展した数値モデルである．数値流体力 学モデルとしてのLBMは，マクロな連続体を仮定したN- S式を離散的に解くのではなく，流体中を並進・衝突を繰 り返す仮想粒子の分布に関する格子ボルツマン方程式を 離散化する．粒子法の一手法として分類されるが，流体 粒子一つひとつの運動をミクロに解かないという意味で，
メゾスケールの数値モデルとして位置づけられる．流体 内に仮想的に配した格子上を運動する粒子の粒子分布関 数を用いて，質量・運動量保存則を満たすように，格子 内の水の密度，水位，流速等の巨視的変数を求める．局 所平衡分布関数（後述）の関数形に応じて，浅水理論式 に対応した長波近似の解やN-S式に対応した解を得ること ができる．
以下では，Zhou（2004）が定式化した浅水理論と等価
なLBMについて概説する．詳細についてはZhou（2004）
や大家ら（2008）を参照されたい．
（2）格子配置
LBMでは空間の離散化のために構造・非構造のいくつ かの格子形状が提案されているが，平面2次元の流れ場に おいて構造格子を用いるのであれば図-1に示す2次元9速 度格子を用いるのが一般的である（例えば，Zhou，2004）． 格子上の粒子の速度ベクトルは式（1）で定義する．
格子ボルツマン法による津波解析コード構築に向けた 実地形適用に関する研究
Applicability of the Lattice Boltzmann Method for Tsunami Modeling
越村俊一
１・村上和幸
２Shunichi KOSHIMURA and Kazuyuki MURAKAMI
The Lattice Boltzmann Method (LBM) is one of the robust CFD models to solve the fluid dynamics. In the present study, the LBM for the Shallow Water Equations is developed to simulate tsunami propagation and coastal inundation.
The offshore boundary conditions and nesting technique for multi-scale grid are developed for practical use of LBM.
The model is validated, for a potential tsunami scenario, by the comparison with the Finite Difference Method. The results imply that LBM has a sufficient accuracy and robustness for tsunami propagation and inundation modeling.
1. はじめに
格 子 ボ ル ツ マ ン 法 （Lattice Boltzmann Method， 以 下 LBM）とは，粒子法として分類されるメゾスケールのCFD モデルである（McNamara・Zanetti,1988 ; Qianら, 1992 ;
Chen・Doolen，1998 ）．流体を，格子上を移動する仮想
な水粒子として近似し，その並進・衝突の時間発展を格 子ボルツマン方程式に基づいて計算し，マクロな流れ場 の諸量（水位，流速等）を求める．LBMでは格子ボルツ マン方程式を粒子の各速度成分の頻度（粒子分布関数f）
を変数として完全に陽的なスキームで離散化し，その解 はナヴィエ・ストークス式（N-S式）と2次精度で一致す ることが数学的に保証されている（渡辺，2006a，2006b）． LBMは単相流だけでなく，気液・液液の混相流の計算に も適した手法であるが，高精度・高分解能の流体計算を 低コストで実現できる点において期待が大きい．
近年，LBMを水害や津波等の災害シミュレーションに 応用しようとする研究が進められており，特にZhou
（2004）が浅水理論と等価なLBM解を提案して以来，潮流 や長波の遡上問題においていくつか適用例が報告される ようになった（石川ら，2004 ; 大家ら，2008）．しかし LBMの実用化に向けては，無反射境界条件や領域接続な どの問題を解決する必要がある．本研究は，LBMに基づ く津波遡上解析コードの構築を目指し，実地形への適用 に必要とされる課題を解決することを目的とする．
本稿では，津波伝播・遡上解析コードとしてのLBMの 実用化と実地形への適用を進める上で，次の3つの課題を 設定し，計算の安定性と精度に関連した検証を行う．す なわち，（1）津波の外洋伝播から陸上遡上に至るまでの
１ 正会員 博（工） 東北大学准教授 大学院工学研究科 ２ 非会員 株式会社NTTデータ
すなわち，全水深h(x,t)は位置xにおける各方向の粒子 分布関数の和で，流速u(x,t)は粒子分布関数の重みに並進 速度ベクトルeαを乗じることで計算される．
（4）局所平衡分布関数
局所平衡分布関数とは，仮想流体粒子の運動が有限な 空間において平衡状態に達した際の粒子分布である．
LBMでは，並進と衝突を繰り返す粒子の各速度方向の頻 度分布が，気体分子運動論における局所平衡解である Maxwell-Boltzmann分布に一致することを局所平衡と定義 し，その時の粒子分布関数を局所平衡分布関数fαeqとして 定式化している．局所平衡分布関数は，時空間（x ,t）に おけるマクロな諸量u= (ui,uj),hを用いて式（7）で求めら れる．
（7）
3. 格子ボルツマン法の実地形適用性
格子ボルツマン法の実用化を考えた場合，i）従来整備 されてきた海底地形データが利用できること，ii）求めら れる精度と計算容量に応じて多様な分解能を選択できる こと，iii）有限な計算領域を扱うための簡便な沖側無反射 透過境界条件を整備することの3点が要件として挙げられ る．海底地形データの利用可能性については，大家ら
（2008）が一様勾配斜面の津波遡上問題において，同一の 空間分解能においてはLeap-frog差分法とLBM解が一致す ることを確認している．以下では，複数の分解能の異な る計算領域をネスティングする多領域接続スキームと沖 側無反射透過境界条件についての検討を行う．
（1）多領域接続スキーム
津波の数値計算には，水深数1,000mの発生域から陸域 まで広大な領域を計算対象とする必要がある．浅海域か ら陸域にかけて高い空間分解能が必要とされる一方で，沖 域（線形長波近似が成立する領域）ではさほど高い分解 能を必要としない場合が多い．従って,沖から浅海域へ空 間分解能を段階的に高めていくネスティング手法を採用 するのが一般的である．いま，図-2に示す格子配置にお いて，大領域の格子点A，Bから小領域の格子点a-dへ1：
3の格子間隔比のネスティングを考える．接続する変数は，
格子点上の全水深h（または水位）と流速uである．大領 域から小領域へ流速uを次式の線形内挿で与え，それを粒 子分布関数に変換する．小文字の添字は小領域，大文字 の添字は大領域の諸量を示す．
…（1）
ここで，eは計算に用いる時間・空間格子を用いてe= Δx/Δtと定義し，eαはα方向の並進速度ベクトル，fαは粒 子分布関数でありeαを持つ粒子の総数を表す．
（3）格子ボルツマン方程式
格子ボルツマン方程式は，衝突演算項に格子BGKモデ ル（Qianら，1992）を用いて，式（2）で表される．以下 では粒子の速度ベクトルの方向を添字αで表し，巨視的な 流れ場を表す変数についてはx，y方向に対してそれぞれ i，jの添字を付す．また次式以降は総和規約を適用する．
…………（2）
ここでfα(x,t)は粒子分布関数であり，時間tにおいて位
置xの格子点上に存在する，α方向の速度ベクトルをもつ 粒子数である．fαeqは局所平衡分布関数（後述），Fは外力 項である．LBMでは，地形変化による効果( zb/ xi)と底面 摩擦項を外力項として式（3）で表現する（Zhou，2004）．
……（3）
ここでnbはManningの粗度係数である．一方，式（2）
のτは単一時間緩和係数で，動粘性係数υを用いて式（4）
で表される．
………（4）
以上の支配方程式から陽的に未知数fαを求め，巨視的 変数である全水深h，流速uを求める．h，uは粒子分布関 数fαの各方向の総和をとることにより，式（5），（6）から 求める．
………（5）
………（6）
図-1 LBM格子の格子配置（2次元9速度格子）
………（8）
全水深hの接続は，式（9）の通り，各格子の静水深H0
を差し引いて水位のみを接続する．
……（9）
上式により，各格子点における諸量から局所平衡分布関 数を求め，それを粒子分布関数として計算を進めていく．
………（10）
一方，小領域から大領域への受け渡しは，空間的に一 致する格子上の水位と流速を直接与えて次式のように求 める．
………（11）
………（12）
最後に，この巨視的変数を用いて局所平衡分布関数か ら各格子の粒子分布関数を求める．
………（13）
（2）無反射透過境界条件
LBMでは，2次元9速度格子の場合，いずれの格子にお
いても9方向の粒子の分布を計算しなければならない．い ま，計算領域内から右側の境界外縁の格子を考える．そ の場合，領域外に出る方向の粒子（図-1の1，2，8番方向）
は並進過程の計算で自動的に排出されるので問題とはな らない．しかし，逆方向（図-1の4，5，6番方向）につい ては粒子の計算を行う必要があるため，右に向かう波条 件下での粒子分布関数を与える必要がある．津波解析を
考えた場合の沖側境界は，水深hが十分に深く，線形長波 近似が可能である．そのとき，式（14）のように進行波 の水位ηと流速uの関係式が成立する．
………（14）
この関係を利用して，4，5，6番方向の粒子分布関数を 次式で求める．
…………（15）
fαeq(h(x,t),u)は同格子点において，全水深h(x,t)，速度 uとしたときの局所平衡分布関数である．uroutは，xにお ける静水深をH0(x)とおいたときに次式から求める．
……（16）
同様に，左側境界において，逆方向の粒子分布関数
（図-1の1，2，8番方向）について次式で求める．また，上 端および下端についての無反射透過条件も同様にして求 めることができる．
………（17）
……（18）
（3）検証
上記手法を一定水深（3,000m），進行波の条件で検証し た結果の一例を図-3に示す．これは，初期水位5mの水面 上昇を図の左下隅（50×50格子分）に設定して計算を実 施し，図の左下から右上に向かって進行する津波の水位 分布を示したものである．図中の□内部がネスティング 領域（小領域）である．微小な水面の乱れが見られるも のの（ネスティング無しの条件で計算した水位差を調べ た結果，±2cm程度であることが分かった），大領域（Δ
x=150m）から小領域（Δx=50m），小領域から大領域への
接続について実用上問題ない精度で計算を行うことがで きる．沖側境界での反射波の発生も抑えることができ，本 手法で実用計算が可能であることが分かった．
4. 格子ボルツマン法の実地形適用に向けた検証
（1）想定津波と地形データ
実際の海底地形データおよび想定津波波源モデルを用
いてLBMを検証した結果を示す．ここでは内閣府中央防
災会議想定の東南海・南海地震津波の波源モデルを利用 し，Leap-frog差分法による計算結果をベンチマークとす る．図-4に示すのは，想定津波の初期水位分布と，検証 計算の領域（三重県尾鷲湾）および波形出力地点（a-cの 3点）である．計算格子は1350mから50mまで1：3で4領
図-2 1：3のネスティング格子配置
域のネスティングを行い，50m格子で計算した尾鷲湾内と 市街地の3地点で計算波形を出力した．遡上計算は50m格 子（尾鷲湾）のみで行い，遡上先端の境界条件は大家ら
（2008）が提案した簡易式を用いた．計算の諸元を表-1に 示す．
（2）差分法との比較検証
図-5に地震発生から22分後（第1波のピークが尾鷲市街 地に到達した時間）における水位分布の比較を示す．上 がLBM，下がLeap-frog差分法（FDM）による結果である．
両者はほとんど一致しており，LBMが従来の手法と変わ らない精度で計算を実行できる．浸水域についても概ね 一致しており，LBMが高い実用性を有していることが分 かる．
また，尾鷲湾内の3地点で出力した計算波形の比較を図- 6に示す．湾口・湾内（a，b点）における波形は，LBMと FDMはほとんど一致している．一方，陸上氾濫の波形（c 点）については若干の差異が確認された．津波遡上先端 部の境界条件と陸上の抵抗則の違いによるものであると 解釈できる．本検証例ではManningの粗度係数n=0.025を 一律に与えているが，LBMでは底面摩擦項を含む外力項 式（3）の計算精度は離散化の仕方により差が生じること が報告されている（Zhou，2004）．ここで示した計算例に は，浅水理論に対して2次近似を満足するCentered Scheme を採用しているが，より詳細な検討が必要であろう．ま た，市街地内でのより詳細な氾濫計算を実施するために，
例えば合成等価粗度（油屋・今村，2002）に対応する抵 抗則も検討する必要があろう．
5. 結論
本研究で得られた結果を以下に列挙する．
LBMに基づく津波解析コードの構築に向け，実地形で の計算を可能にするための多領域接続計算法と無反射境 界条件を整備し,その妥当性を検証した．従来整備されて きた海底・陸上地形データを有効に利用し，大領域から 小領域，小領域から大領域への接続がスムーズに処理で きることが実証できた．
沖側境界の無反射境界条件を，線形長波の進行波の条 件を用いて領域外から内側に向かう粒子分布関数を疑似 的に与えることで設定し，沖側境界での反射波の発生も 抑えることができた．
LBMの実地形への適用として東南海・南海地震津波を 想定し，三重県尾鷲市の津波遡上解析を実施した．従来 の差分法との比較を通じてLBMの妥当性を検証した結果，
LBMによる津波の遡上解析を実地形においても実用上問 題ない精度で実施できることがわかった．一方，陸上の 津波遡上計算については，遡上先端部の境界条件や底面 摩擦等の抵抗則など詳細に検討する必要がある．特に，建 物や構造物の影響をどのように考慮するかは今後の課題 であろう．
図-3 進行波の条件における領域接続スキームとの無反射透過境界条件の検証結果（□が小領域）
図-4 検証計算の領域と想定東南海・南海地震津波の波源
（実線は隆起，点線は沈降，コンター間隔は0.5m）
計算法 支配方程式 空間格子間隔 時間格子間隔 単一時間緩和係数（τ）
潮位条件 対象地震 境界条件
現象再現時間
LBM / Leap-frog 差分法 格子ボルツマン方程式/浅水理論式
1350 / 450 / 150 / 50m 0.5sec.
0.500-0.503 平均潮位（TP+0.0m）
中央防災会議想定東南海・南海地震 沖側：自由透過
陸側：完全反射および遡上（尾鷲湾のみ）
地震発生から1.5時間 表-1 検証計算の諸条件
ツマン法による浅水長波流れ解析，第18回数値流体力学シ ンポジウム講演概要，D4-4.
大家隆行・越村俊一・荒木 健（2008）：格子ボルツマン法に 基づく津波遡上シミュレーション手法の開発，海岸工学論 文集，第55巻，pp.221-225.
渡辺 正（2006a）：格子ボルツマン法（1），ボルツマン方程 式 か ら 格 子 ボ ル ツ マ ン 方 程 式 へ ， 応 用 数 理 ，Vol. 16，
No.1，pp.31-35.
渡辺 正（2006b）：格子ボルツマン法（2），ボルツマン方程 式からナビエ-ストークス方程式へ，応用数理，Vol. 16，
No.2，pp.64-69.
Chen, S. and G. D.Doolen (1998)：Lattice Boltzmann Method for fluid flows, Annual Review of Fluid Mechanics,1998 Vol.30, pp. 329-364.
McNamara, G. R. and G. Zanetti (1988) : Use of the Boltzmann Equation to Simulate Lattice-Gas Automata, Physical Review Letters, 61, pp.2332-2335.
Qian, Y. H., D. d Humieres and P. Lallemand (1992)：Lattice BGK Models for Navier-Stokes Equation, Europhysics Letters, 17 (6), pp.479-484.
Thürey, N. (2007) : Physically based Animation of Free Surface Flows with the Lattice Boltzmann Method, Universtität of Erlangen-Nürnberg, Ph.D. Thesis, 145p.
Zhou, J. G. (2004)：Lattice Boltzmann Method for Shallow Water Flows, Springer, 112p.
LBMは，本稿で報告した浅水理論と等価なモデルに加 え，自由表面流れの解析も可能な非静水圧3次元計算への 拡張も進んでいる（たとえばThürey，2007；荒木・越村，
2009）．沿岸域から複雑な街路を含む市街地氾濫の統一的
な解析を可能にする2D-3Dハイブリッドコードへの展開 も今後の課題である．
謝辞：本研究の一部は，科学研究費補助金（挑戦的萌芽，
代表：越村俊一，課題番号：21651078）の補助を受けて 実施された．ここに記して謝意を表する．
参 考 文 献
油屋貴子・今村文彦（2002）：合成等価粗度を用いた津波氾濫 シミュレーションの提案，海岸工学論文集，第49巻，
pp.276-280.
荒木 健・越村俊一（2009）：格子ボルツマン法による自由表 面流れの解析，海岸工学論文集，第56巻，pp.56-60.
石川裕士・立石絢也・田中聖三・樫山和男（2004）：格子ボル 図-5 尾鷲湾における東南海・南海地震津波の水位分布（地 震発生から22分後，上：LBMによる計算（本研究），
下：差分法による計算）
図-6 尾鷲湾内と市街地内の計算波形の比較（a：湾口，b：
湾奥，c：市街地内）
