Expanding Slice Grid 機能付き粒子法による広域・高詳細３次元津波遡上解析の性能評価
浅井光輝（九州大・工）西浦 泰介（海洋開発研究機構）古市Mitsuteru ASAI, Kyushu University
Daisuke NISHIURA, JAMSTEC
Mikito FURUICHI, JAMSTEC

幹人（海洋開発研究機構）

The Tsunami caused by the 2011 Tohoku Earthquake caused tremendous damage. We should prepare
carefully for future huge earthquakes such as the Nankai trough earthquake in order to prevent damages fromthese earthquakes and tsunamis. Due to future tsunami disaster prevention planning, it is strongly desired asimulation tool that accurately predicts possible damage level for structures and human resources in theinundation area. To achieve the goal, we proceeded the development of large scale tsunami run up simulation
using 3-dimensional particle method, which is capable for accurate prediction of large scale inundation area.
In this report, a high fidelity tsunami run-up model in the part of Fukushima nuclear power plant, which canresolve the inside of structure of the power plant, is simulated using the expanding slice grid method.
1. はじめに
東日本大震災では観測史上最大の大津波が発生し，東北地方の沿岸部を中心に，甚大かつ広域的な被害をもたらした．最大クラスを想定した災害への備えが不十分であったことが指摘されており，想定し得る最大規模の浸水等へ対策の措置が講じられている．減災のためには，津波による沿岸地域への安全性・危険性を把握することが重要であり，特にソフト防災の一環として行う津波防災地域づくりにおいては，浸水域が的確に推定できる手法が必要とされている．本研究では，可能な限り正確に津波被害を想定するための，特に遡上後の建屋への浸水までを含む津波解析手法を開発することを目的とし，高詳細な粒子解析モデルを使った大規模並列計算を九州大学スパコン Ito の全ノードを使って実施した．
２．動的負荷分散機能付き陽的 ISPH 法
解 析 手 法 と し て 安 定 化 ISPH(Incompressible Smoothed
Particle Hydrodynamics)法 を陽解法化した陽的 ISPH 法を用いることでメモリと計算コストの両者を削減することで限られた計算資源にて大規模計算を可能とし，また拡張ス
ライスグリッド法による動的負荷分散技術と組み合わせることで高並列計算環境においても効率的な計算コードを開発した．
２．１ 安定化 EISPH 法Asai1)らは，非圧縮性流体解析用の ISPH 法において，解析精度と安定性の両者の向上を目的とし，微圧縮性を一部許容した下記の圧力ポアソン方程式を提案した．
∇2 𝑝𝑖𝑛+1 ≈∇2 𝑝𝑖𝑛+1 =𝑗𝑚𝑗

𝜌0𝜌𝑖0 − 𝜌𝑖𝑛∆ ∙ 𝑢𝑖∗ + 𝛼∆𝑡∆𝑡 2(1)

2 𝑟𝑖𝑗 𝛻𝑊 𝑟𝑖𝑗 , ℎ(𝑝𝑖𝑛+1 − 𝑝𝑗𝑛+1 )𝜌0 𝑟𝑖𝑗2 + 𝜂2(2)表記の簡便化のため，Aij，Bi を，圧力ポアソン方程式の係数行列とソース項とすると次式により記述できる．

𝑝𝑖𝑛+1 =

𝐴𝑖𝑗 = 𝑚𝑗

𝐵𝑖 +𝑛+1
𝑗 𝐴𝑖𝑗 𝑝𝑗

(3)

𝑗 𝐴𝑖𝑗

2 𝑟𝑖𝑗 𝛻𝑊(𝑟𝑖𝑗 , ℎ)
𝜌0 𝑟𝑖𝑗2 + 𝜂2(4)

𝐵𝑖 =

𝜌0
𝜌𝑖0 − 𝜌𝑖𝑛∆ ∙ 𝑢𝑖∗ + 𝛼∆𝑡∆𝑡 2(5)式(3)において，n ステップから n+1 ステップまでの微小な時間間隔では圧力の変動が小さいものと仮定し，右辺中の近傍粒子の圧力は前ステップで近似できるものとすれば，圧力が陽的に更新できる．
𝑝𝑖𝑛+1 =𝐵𝑖 +𝑛𝑗 𝐴𝑖𝑗 𝑝𝑗𝑗 𝐴𝑖𝑗(6)

以上，安定化 ISPH 法における圧力の陰的評価箇所を簡略化し陽的な評価を行う手法を安定化 EISPH 法と呼ぶ．
２．３ 2 次元スライスグリッドによる動的負荷分散
上記の安定化 EISPH 法により，従来型の半陰解法と比べると大幅にメモリと計算量を削減でき，大規模な計算が実現可能となった．高並列な計算機環境において，同コードによる並列計算を効率的に実施するには，動的負荷分散を考慮した領域分割型の並列解析機能が必要となる．本研
究の対象は津波氾濫解析
Fig. 1 2D スライスグリッド法であるため，水深方向（鉛直方向）には浅く，広領域な解析が必要となるため，鉛直方向には領域を分割せず，平面方向のみを Fig. 1 に示す 2 次元スライスグリッドより領域分割を行った．領域分割方向を xy 平面とすると，まずは x 方向に(I=1,2,3,…Nx)の列要素に一次元領域分割を行い，その後に y 方向には各列ごとに最適な分割(J=1,2,3,…Ny)を行うのがスライスグリッド法である．各要素は袖領域を持ち，その領域の粒子群は MPI 通信を行う．
粒子法による計算では，時間とともに粒子配置が大きく変動するため，分割される領域を更新しながら動的に最適な負荷分散を行う必要がある．またこの種の粒子法の並列計算では，負荷分散の指標として粒子数を選び，粒子数を均等化することで負荷を均等化することが多い．しかし，自由表面近傍，壁粒子近傍，内部領域では計算コストが異なり，また通信量も無視できない．そこで，本研究では通信を含んだ計算時間を均等化することにし，これを偽ニュートン法 による動的負荷分散手法により実現した．

２．４ 拡張スライスグリッド法
粒子法では，固体境界を境界粒子としてモデルすることが多い．津波遡上現象のような問題を扱うには，時間とともに解析領域が広がり，これに応じて必要となる壁粒子の領域も拡張していく．またこの浸水領域は，一般的には解析前には事前に特定できないため，遡上の可能性がある領域すべての壁粒子を記憶しておく必要がある．これをすべてメモリ上に記憶しておくことは資源の無駄であり，有効となる壁粒子のみを選択しながら，計算が進行するにつれて必要な壁粒子データを更新していくことが大規模計算では必須である．そこで本研究では，前述の 2 次元スライスグリッドは水を含む領域のみを対象とし，津波遡上に併せて必要となる領域を拡張しながら，分散型の並列計算環境に事前に均等に配分しておいた壁粒子の情報を逐次読み込みながら解析をする拡張スライスグリッド法を提案し，その有効性を計測することにした．
３．解析対象本研究では建物内部への浸水解析を行う事を目的に，2011 年東日本大震災による津波で浸水した福島第一原
子力発電所を解析対象として設定し，.
1 号機建屋内部の浸Fig. 2 福島原発第 1 建屋モデル
水解析を行った．対象モデルの大きさは Fig. 2 に示すように 326m × 537m と設定した．
建屋内部の構造までを詳細にモデル化するため，粒子間隔は 30 ㎝が必要最低限となる．
建屋内部の詳細図はFig. 3 に示す．このときの粒子数及び解析条件をまとめてTable.1 に示す．また，EISPH 法を採用したため，計 5000万粒子（津波の流入を始めると 1 億以上の粒子）の問題が約 50 ノード程度の計算機環境で解析できるコstairードである．この先端的計算科学研究プロジ
ェクトの実施前に，京都大学のスパコン・システム A(Camphor2)を使い，50 ノードの計算inflow機(100MPI プロセス)で計算実行した．このプロジェクトではさらにFig. 3 建屋内部構造大規模な問題を設定するため，粒子間隔を 10 ㎝にした詳細な解析モデルを作成した．これは初期で総粒子数が 4 億以上，浸水まで解析すると 10 億近くの粒子数となる．
このモデルを九州大学スパコ
ン Ito システム A の全 2000 ノードを使用し，半数の 1000ノード使用時の場合と計算時間を比較検討した．

particle

Table. 1

30 ㎝粒子間隔モデル

wall

23,710,256

water

28,500,186

total

53,088,027

dt

0.0005α0.0002

４．解析結果10cm モデルの解析結果は可視化まで実施できていないため，ここでは参考までに 30 ㎝モデルの解析結果を Fig.4に示す．Fig.4 a)には浸水の発展の様子と各時間ステップでの領域分割を示している．本解析モデルでは，遡上現象だけでなく建屋の内部への浸水までを含んだ高詳細なモデルであり，建屋へ浸水した水は地下まで浸水する複雑な現象を解いている（Fig.4 b)を参照）．安定化 EISPH 法により計算負荷を低減したアルゴリズムを用い，さらには拡張スライスグリッドによる動的負荷分散技術を併用することで，50 ノード程度の計算資源での解析を可能とした．
2019 年度前期の先端的計算科学研究プロジェクトでは，10cm 解析モデル(初期粒子数は 4 億以上)を 2000 プロセスに分割し，2000 ノードで無事に計算できることまでは確認できたが，計算モデルの修正，解析パラメータのチューニングに時間を要してしまった．このため，2019 年度後期の同プロジェクト課題を継続申請し，計算効率を計測することにした．

Time：27.25sa)b)

Time：158.5s

津波氾濫解析結果

建屋内部への浸水の様子

Fig. 4 解析結果例(参考：30 ㎝モデル)

５．参考文献1) 小笠原圭太，浅井光輝，古市幹人，西浦泰介： 大規模津波遡上解析に向けた安定化 ISPH 法の 陽的時間積分ス
キーム，土木学会論文集 A2， Vol 73(2), I_397-I_404, 20172) M. Asai,A.M. Aly ,Y. Sonoda , and Y. Sakai: A stabilized
Incompressible SPH method by relaxing the Density invariancecondition, Journal of Applied Mathematics, vol.2012,139583,pp1-24, 2012