学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

14-NA25

沈み込み帯の巨大地震を対象とした大規模並列地震波・津波伝播シミュレーション
竹中博士（岡山大学）
概要

本課題研究では時間領域差分法にもとづいて、沈み込み帯の巨大地震への応用を目的とした大規模地震波シミュレーション手法に関する研究に取り組んだ。
今年度の主な研究内容は次の通りである。［1］GPU演算カーネルの性能向上（ 5.1）：TSUBAME-2.5で実行するプログラムのGPU演算部分についてNVIDIAのKeplerアーキテクチャに対応した最適化を検討した。［2］東北地方太平洋沖地震の解析に向けた短周期地震波再現性の検討（ 5.2）：震源自体の複雑さを回避できる小地震データを対象として、TSUBAME-2.5を用いた大規模地震波シミュレー
ションによってグリーンテンソル波形を生成した。それらを用いた逆問題解析によって、特に短周期の表面波の再現性について検討を行った。［3］摂動カーネル計
算の効率化（ 5.3）：波形再現性を向上させるためには構造モデルの改良が必要である。そこで構造モデル摂動に関する波形の感度を表す摂動カーネル計算のためにプログラムを改訂し、TSUBAME-2.5において演算性能に大きな影響を及ぼすことなく出力ファイルサイズを大幅に低減することができた。［4］東京大学FX10での計算高速化（ 5.4）：各種の最適化を施したプログラムを利用して、FX10での実行性能が改善されることを確認した。この結果、これまではFX10のshortジョブクラスでは実行不可能だった大規模計算が実行可能となった。

１．共同研究に関する情報(1) 共同研究を実施した拠点名
国立大学法人 東京大学
（Fujitsu PRIMEHPC FX10）

氏名

所属

役割

青木尊之

東京工業大学・

副代表・GPU計算の最学術国際情 報セ 適化、大規模並列化、

岡元太郎

国立大学法人 東京工業大学（TSUBAME 2.5）

中村武史

(2) 共同研究分野
小松正直■ 超大規模数値計算系応用分野□ 超大規模データ処理系応用分野□ 超大容量ネットワーク技術分野

ToshiroTanimoto□ 超大規模情報システム関連研究分野
(3)

ンター

計算結果の可視化等

東京工業大学・副代表者・プログラム

理工学研究科

開発と応用計算海洋研究開発機共同研究者・CPU版プ構

ログラムの開発と応用

岡山大学・大学

共同研究者・構造モデ

院自然科学研究

ル構築と地震波・津波

科

計算

カリフォルニア

共同研究者・短波長波

大学サンタバー

動場の物理に関する基バラ校

礎的考察・摂動カーネルの計算と考察

参加研究者の役割分担
氏名

所属

役割

竹中博士

岡山大学・大学

研究代表者• プログラ

院自然科学研究

ム開発と応用計算、全科

テーマの総括２．研究の目的と意義
2011年3月11日に発生した東北地方太平洋沖地震（マグニチュード９）は、強い地震動と巨大な津波によって東日本地域に計り知れないほどの地震災害・津波災害をもたらした。この地域ではマグニチュード９に達するほどの巨大地震が発生することはほとんど想定されていなかっ!1

学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

たため、この地震の詳細な断層破壊過程を明らかにし、巨大地震発生に至った過程や強震動・津波の生成メカニズムを考察することが重要な研究課題となっている。このことは他の地域の海溝型巨大地震による強震動や津波を評価・予測するうえでの基盤ともなるので、防災面でも重要な課題であると位置づけられる。

とづく大規模地震波シミュレーション手法と地震
波・津波統合シミュレーション手法を利用して、海溝型巨大地震の研究に関する次の項目に取り組んだ。(1) 3次元構造モデルの検証の高度化［平成25年度最終報告書 ３.２］：東工大TSUBAME-2.0を用いた大規模並列GPU計算により摂動カーネルを求めて可視化を実施した。それらの結果に基づき3次元不均質構造モデル内での短周期地震波の伝播の物理に関して考察を進め、論文を投稿した。(2) 東北地方太平洋沖地震の破壊過程［同 ３.３］：破壊過程の逆問題解析の分解能や信頼度を検証するために、逆問題解析に関するシミュレーションを実施
し、グリーンテンソル波形の違いが解析結果にどのように影響するかについて考察を行った。波形計算に用いる構造モデルを改良し、実データにもとづく破壊過程解析も進めた。(3) 地震波・津波統合シミュレーション［同 3.4］：南西諸島地域の3次元構造モデルを用いて、東大FX10を用いた地震波・津我々はこのような課題に取り組むために大規模
地震波・津波伝播シミュレーション手法（時間領域差分法FDTD: Finite-Difference Time Domain）の開発と応用を進めてきた。この方法では食い違い格子を用いて媒質を離散化し、格子点に粒子
速度ベクトルの各成分（合計3成分）と応力テンソルの独立な各成分（合計6成分）とを変数として割り当てる（図１）。

平成25年度の課題では東北沖地震の破壊過程の推定や地震波・津波統合シミュレーショ
ンなどを主な目的として研究を行った。平成26
年度にはこれらの研究を進展させるとともに、演算能力やメモリが増強されたTSUBAME-2.5を利用してこれまでよりも短波長（短周期）の現象の再現性について研究を進める。

波統合プログラムによる地震動計算を行った。

５．今年度の研究成果の詳細
以下の項目のうち 5.3、 5.4が中間報告以降の成果である。

5.1 GPU演算カーネルの性能向上

図1 食い違い格子の単位

TSUBAME-2.5を用いる地震波シミュレーショ
ンプログラムのうち GPU演算部分（GPUカー

セル[Okamoto et al.
2010]。Viは速度成分、

ネル）についてNVIDIAのKeplerアーキテクチャ（K20X）に対応した最適化を検討する。
なお検討に先立ってTSUBAME-2.5の開発
環境が更新されたので、その効果を確認しておく。昨年度の研究で最適化以前の約1.5倍の性能を達成した［岡元, HPCI成果報告書,
2014］。その時点でのTSUBAME-2.5の開発
環境はCUDA-5.0であった。その後、開発環境
がCUDA-6.0に更新されたので（2014年8月）、
全く同じソースコードを用いて性能を比較する。図2に示すようにCUDA-6.0への更新により本プログラムでも性能が約18％向上し、1GPUで216.1GFlopsに達することになった。

τijは応力テンソル成分を示す。

３．当拠点公募型共同研究として実施した意義本研究では多数のGPUを用いた大規模並列計
算によって地震波伝播のシミュレーションを行う 。 東 京 工 業 大 学 国 際 学 術 情 報 セ ンタ ー のTSUBAMEスーパーコンピューターは2013年にGPUを更新して倍精度5.7ペタフロップスのピー
ク性能を有する、世界でもトップクラスのGPUクラスタである。そのため、本研究を進めていくうえで最適な環境を利用することができる。
また、京コンピュータ等のCPUクラスタシステ
ムの環境への対応も必要であることから、東京大学のFX10ではCPU版プログラムにより地震波•津波統合シミュレーションを実施する。

700

GFLOPS

600
500
400
300
200

４．前年度までに得られた研究成果の概要

100

前年度の 課題では本年度と同じく、差分法にも

0!2

CUDA-5.0
CUDA-6.0

1GPU

2GPU

3GPU

学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

図2（前ページ） 地震波シミュレーションプロ

16 14
16 28
16 42
16 56

240

グラムのベンチマーク。同じソースコードを
CUDA-5.0とCUDA-6.0でコンパイルして性能を比較したもの（副領域サイズ 320 320 160、200

GFLOPS

格子間隔 0.2km、時間間隔 0.01s、ステップ数
2000）。共有メモリは48kB、スレッドあたり
レジスタ数は64個、物性定数配列に共有メモリを適用した。非弾性要素数は6個なので、1単位セルあたりの場の量の変数は45個となる。ここ

160

120

の性能は、差分法演算ループを構成するGPUカーネルの繰り返しに要した総時間と、

80

CUDA-6.0のプロファイラ（nvprof）で測定した演算回数（単精度）をもとに計算したもので

32

64

96

128

160

レジスタ数

ある。したがって構造モデル読み込み・設定な
どの主要な演算前後の処理時間は含まない。ま図3 1個のGPUでブロック数とスレッド数とをた演算中のファイル出力時間も含まない。測定

パラメータとしたベンチマーク結果。ブロック

にはTSUBAME-2.5の1ノードを使用した。

数の違いにより領域サイズも図2の場合とは異なる値になる。それら以外の設定は図2との場

5.1.1 ブロック個数とレジスタ数の見直し
我々が開発したGPU版の地震波シミュレー
ションプログラムは3次元媒質を取り扱うが、GPUカーネルでは2次元ブロックを採用している。これはオンチップの共有メモリが小さく3次元ブロックを収めるのが難しいためである［Okamoto et al. 2010, 2013］。この2次元ブロックの個数は、図2に示したこれまでの計
算では400個としていた（32 8のブロックが10 40個）。しかしTesla K20X に使われているGPU（GK110）では、16個のブロックを同時に割り当てられるマルチプロセッサ（SMX）
を14個備えているので、ブロックの総数を
16 14＝224個の倍数にした方が性能が上がる可能性がある［河合、2013］。そこで、前項のプログラムのブロック個数を224の倍数に変更することを試みる。すなわち32 8のスレッ
ドブロックサイズを固定して（1スレッドが1単位セルに割り当てられる）、このブロックの個数が16 14（224個）、16 28（448個）、
16 42（672個）、16 56（896個）の4通りの場合の性能を検討する。
ここでは同時にスレッドあたりのレジスタ数も検討する。レジスタが多いほどスレッドの性能は上がるが、総レジスタ数には上限があり、レジスタ数と、一度に処理できるスレッ
ド数とにはトレードオフがあるためである。
K20Xではスレッドあたりの最大レジスタ数が225個までに引き上げられたので、32、64、96、128、160個のレジスタ数を検討する。

合と同じである。演算回数にはプロファイラによる測定値を用いた。CUDA-6.0を用いた。

図2に示したベンチマーク結果によると、本研究のプログラムでもレジスタ数の選択は性能に大きな影響を及ぼすことがわかった。一方でブロック数の違いによる性能への影響は比較的に小さいこともわかった。
なお、上述の条件のもとで最大性能を与えるのはブロック個数が16 14（224個）、レ
ジスタ数が64個の場合で、226.6 GFLOPSである。これはブロック個数を最適化していな
い図2の場合の性能（216.1 GFLOPS）から約5％の向上である。
5.1.2 Read-Only キャッシュ
KeplerではL1キャッシュメモリにキャッシュ
されるのはローカル変数のみであり、グロー
バルメモリの変数はL2キャッシュに入る仕様になった。その代わりに書き換え不要の変数を比較的に簡単な手続きでキャッシュできる
Read-Onlyキャッシュ（ROキャッシュ）が用意された。本研究で使うプログラムでも物性
定数配列は演算中に書き換えないので、これらをROキャッシュに割り当てることを検討する。前項までの結果は物性定数配列を共有メ
モリに割り当ててコピーする方法を使っているので、共有メモリ、ROキャッシュ、グロー
バルメモリ変数をそのまま使用（L2キャッシュに入る）の3通りを比較する。テストには前項!3

学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

での最大性能ケースのパラメータを用いる。

性能

SM

GM

RO

RO＊

226.6

220.2

211.7

87.1

震動波形を同時に用いた逆問題解析により、最適震源位置、モーメントテンソル、震源時間関数を推定して、波形の再現性を考察した［Okamoto et al., AGU 2014 他］。ただしこの考察には以下のような制約条件があった。

表1 1個のGPUの場合に物性定数配列の操作を① 震源位置の探索をプレート境界に限定② やや長周期側（12-40秒）の波形を利用変えたときのベンチマーク結果（性能はGFLOPS）。SM：共有メモリを使用。GM：グロ ーバルメ モ リ を 使 用 （ L 2 キ ャ ッ シュ に 入る）。RO：ROキャッシュを利用。RO＊：RO
キャッシュを使うが共有メモリサイズも16kBに変更。他の設定は 5.1.1の最大性能の場合（14 16、レジスタ数64)と同じである。
CUDA-6.0を用いた。

表1からわかるように、本研究のGPUカーネ
ルではROキャッシュを使うとL2キャッシュの場合よりも性能が低下する。ただしこれは今回のベンチマークで使った構造モデルが単純で物性定数配列サイズが数要素と小さいため、
L2キャッシュが有効に働いた可能性も考えられる。一方で物性定数について共有メモリを利用しない場合でも極端な性能低下は見られなかったことから、複雑な構造モデルで物性定数が多くの要素を必要とする場合でも共有メモリを消費せずに済ませられる見通しも得られたと言える。

図4 東北地方太平洋沖地震の解析で用いる強震動観測点（▲）と、本報告で用いる小地震（2003-Mw5.8、小さな★；2003年11月1日22時10分、Mw5.8）を示した図。なお2003年の小 地 震 の モ ーメ ン ト テ ン ソ ル 解 （ G l o b a l

5.1.3 共有メモリとL1キャッシュのサイズ
本研究のGPUカーネルではローカル変数を多用する。そこでL1キャッシュを大きくする
ことで性能が向上するかどうかを見るために、
上記のROキャッシュの場合に共有メモリを48kBから16kBに縮小し、L1キャッシュを16kBから48kBに拡大した例を検討した。表1に見るように、この場合には性能が著しく低下することが確認され、本プログラムでは共有メ
モリが性能に重要であることを確認できた。

CMT）・2011年東北地方太平洋沖地震の震央（大きな★）・同地震のGlobal CMT解（GCMT）とW-phase CMT解（WCMT）・同
地震の解析で用いる断層モデルを構成するノー
ド点（●）もそれぞれ図示した。

上記の①については、この小地震がプレート境界型の発震機構解を持つことから探索対象となる位置を限定する必要がある場合には妥当な条件である。しかし仮定した構造モデルは完全ではないため、①の条件を外したときの最適位置から構造モデルを改良するための情報が得られると考えられる。②について
ここでは、やや短周期側（7-40秒）の波形も利用して検討する。わずかな周期帯の違いのようであるが、海域に特徴的な短周期表面波が顕著になる周期がおよそ10秒程度以下であることから、後に見るように波形には大きな違いが現れる。

5.2 東北沖地震の破壊過程の検討
5.2.1 小地震による波形再現性の検討の背景

平成25年度に実施した破壊過程の推定結果を用いて、東北沖地震の短周期地震波の再現性を検討する。ここではまず、震源の複雑さを回避できる小地震の観測波形の再現性を検討する。
我々はこれまでにこの小地震（2003年
mw5.8）について、遠地実体波波形と近地強!4

学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

5.2.2 大規模地震波シミュレーション解析に使うグリーンテンソル波形を生成するために、大規模地震波シミュレーションを実施した。これは本項目の解析では、 5.2.1に述べたように、プレート境界から外れた点に対応するグリーンテンソル波形が必要となるためである。今回は暫定的に下記のパラメー
タで2観測点（IWTH04とMYGH06、図4）の水平動速度成分を計算した。合計4成分なので、TSUBAME-2.5を利用して4回の大規模シミュ
レーションを実施した。
設定項目利用したパラメータ

全格子サイズ

5120 3200 1440

副領域サイズ320 320 480

格子間隔0.15 km

全領域サイズ768 480 216 km

時間間隔

0.0071429 s

非弾性要素数

3

GPU数

480

ステップ数

38094 3、23338 1

総実行時間

8777.1 s

深さと一致する。
これらの解析で得られた波形を図5に示す。
長周期側（図5左側）では波形全体にわたって
理論波形が観測波形をおおむね再現できていることがわかる。一方、短周期側（図5右側）
では、50秒付近までの主要な波はある程度再現できているものの、後半に到来する短周期の波を再現できていないことがわかる。これは、海水層や厚い堆積物層の存在によって発生する短周期の表面波であると考えられる。
これらの結果、つまり短周期表面波が励起されていないことと最適震源深さはほぼ構造モデルのプレート境界にあること、および海の水深は比較的によくわかっていること、などを考慮すると、構造モデルのうちの堆積物
構造を改良していく必要があるという見通しが得られる。

表2 グリーンテンソル波形計算のための大規模
地震波シミュレーションパラメータ。総実行時
間はステップ数が38094である3回の平均値を記載した。
図5 格子探索で得られた最適震源による理論

5.2.2 格子探索と小地震の波形再現性

波形（SYN）と観測波形（DATA）の比較。

最適震源位置を推定するために、小地震の気象庁震央位置（図4の★）を囲むように、水平方向各4点（約10km間隔）、深さ方向に10点（0.6km間隔）、合計160個の候補震源位置からなる格子を設定し、それぞれの震源位置
で暫定的な非線形逆問題解析を行った［Okamoto & Takenaka 2009 の方法］。
その結果、長周期側（12-40秒）では気象庁
震央位置よりも北北東に約10kmずれた位置が最適位置として得られた。深さは10.3kmであり、構造モデルのプレート境界（11.5km）よりも若干浅い。一方、短周期側（7-40秒）では気象庁震央よりも少し（約6km）陸寄りで、深さ11.5kmとなった。これはプレート境界の

左側：長周期波形による結果。右側：短周期波形による結果。

このような直観的な見通しを定量的に検証する方法が、次節に述べる摂動カーネルの方法である。

5.3 摂動カーネル計算の効率化
上述のように、地震波波形の再現性を検討するには、構造モデルの波形への影響を定量的に調べる方法が必要となる。そのための一つの方法が、昨年度（平成25年度）にも研究を行った、構造摂動に対する波形の感度を積!5学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

分カーネルとして求める方法である（以下では 「 摂 動 カ ーネ ル 」 と 呼 ぶ ） ［ T a r a n t o l a
1984; Tanimoto 1990; Tromp et al. 2005 など］。手法の概要は昨年度（平成25年度）の最終報告書に記載した。なお申請書にも記載したように、この項目では平成25年度に引き続いて、人工データが充実しているカリフォル
ニア州を対象としてカリフォルニア大学の研究者とも協力して研究を進めた。
この計算では計算領域内の各点での速度・応力成分を出力する必要があるため、出力デー
タが膨大な大きさになる。そのため昨年度は1枚の断面だけに限定して出力するようにしていた。それでも計算例では出力ファイルサイ
ズは 1.5 TB ほどに達した。このことは摂動
カーネルの研究を進めるうえで大きな制約になっていた。
ここで、出力ファイルサイズが大きくなるのは時刻歴波形を全て保存していたためである。しかし摂動カーネル計算で必要になるのは時刻歴波形全体ではなく、特定の周波数に関するそのフーリエ成分である。
設定項目

利用したパラメータ

全格子サイズ

4160 2560 1000

副領域サイズ320 320 200

格子間隔

0.08 km

全領域サイズ333 205 80 km

時間間隔0.004 s

非弾性要素数

3

GPU数

520

ステップ数

40000 2例

出力格子点416 256 200

出力周波数の個数

10

出力ファイルサイズ

20 GB

総実行時間

9867.8 s, 10096.9 s

周波数成分のみを出力するような改訂を行った。この改訂により出力ファイルサイズを大幅に低減できた。以下に紹介する例では、1周波数出力の場合にはもとの時刻歴波形全デー
タの数千分の1の出力サイズとなる（表3）。
この計算の一例として、南カリフォルニアを対象とした大規模計算例を紹介する。これは昨年度のJHPCN報告書および Tanimoto &Okamoto (2014) の計算条件を用いた例である。本計算は東工大のTSUBAME-2.5を用いて実施した。計算パラメータは表3に記載した。
2回の計算が記載されているのは、摂動カーネ
ルの計算において人工起震点（図6のSRC）および観測点（図6のREC）のそれぞれに置いた
震源からの波動場を干渉計算させるためである。図6には摂動カーネルの例を示す。ここでは深さスライスを10kmおきに3枚示した。しかし計算では深さ方向に400mおきに出力したので、より細かな検討も可能である。

表3 摂動カーネル（図6）を計算するための差

図6 構造モデルの摂動に対する波形の感度を示

分法パラメータ。

した積分カーネル（摂動カーネル）。対象周波数1.1Hzに関する剛性率の摂動カーネル実部を示した。振幅が大きい位置ほど波形に大きな影

そこで今年度は、出力ファイルサイズを圧縮する方法をプログラムに実装することを進めた。具体的には、差分法演算中は時刻歴波
形を配列に保存しておき、差分法演算ループ
終了後に波形にFFTを施して、最終的に特定の響がある。 （SRC）は人工起震点を置いた位
置、 （REC）は観測点位置である。一番上は地形図である。摂動カーネルの振幅は規格化して図示した（ただし一部は飽和している）。

!6

学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

図から、浅い場所では起震点と観測点を直接に結ぶ部分の感度が高いこと、その部分の幅は1.1Hzという短周期でも20km程度に達す
ること、起震点の後ろ側（ロサンゼルス盆地など）では波形の振幅が大きくなることに起因して摂動カーネルもかなり大きくなること、深くなると摂動に対する感度が次第に弱くな
ること、などが定量的に把握できる。
このような摂動カーネルを系統的に得ることができれば、構造モデルを改良することも可能となる。
図7 淡路島地震の震央（☆印）と、サンプル波

5.4 東京大学FX10での計算の高速化
これまでに開発してきた地震波・津波統合
シミュレーションプログラム［Takenaka et
al. 2012］では、津波まで含めると典型的な時間ステップ数が10万を超える。そのため計算時間が長くかかり、東京大学FX10で大規模ジョブを実行しづらい（待ち時間が長くなる）という実用面での問題があった。そこで今年
度はFX10においてプログラムの実行速度を高速について検討した。プログラムに導入された最適化の主な項目は次の通りである。
• 非同期通信による通信時間の隠蔽• 時間発展ループ内での分岐やサブルーチ
ン呼び出しの削減• 配列融合によるメモリアクセス削減• ループアンロールとパイプライン化
なおこれらの最適化は、これまでのところ地
震波のみの計算の場合に施されている。津波
項を含む場合についての最適化は平成27年度に継続して実施することとした。
設定項目利用したパラメータ

全格子サイズ

1601 2401 401

格子間隔

水平 0.558 km、鉛直 0.5 km

全領域サイズ892.8 1339.2 200 km

時間間隔

0.025 s

非弾性要素数3

ノード数

70 （1120コア）

ステップ数

20000

総実行時間

3.18 時間形を出力した観測点（△印）。

図8 淡路島地震のシミュレーション結果例。図7の観測点における地動の上下動速度成分を震央距離に従って並べて描いたもの。安定に計算できることがわかる。

これらの最適化によって高速化されたため（計算条件にも依存するが、津波項の計算を含む旧プログラムと比較した例では約20倍程
度）、これまで東京大学FX10のshortジョブク
ラスでは実行できなかった規模の計算が、同ジョブクラスで実行できるようになった。ここでは、そのような大規模計算の一例を紹介
する（表4）。この例は上述のように地震波のみの計算に限定する（津波項を含まない）。
また既存結果との比較計算のための例である
ことから、沈み込み帯地震ではないが、2013
年淡路島地震を対象とした（図7）。計算には3次元構造モデルを採用し、陸海地形や海水層も含む。サンプル波形の出力結果を図8に示す。

表4 淡路島地震（2013年4月13日、M6.3、図7）の地震波シミュレーションパラメータ。東京
大学FX10のshortジョブクラスで実行した。

!7

学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

この例のように、本研究では西日本全域にわたるような領域で、しかも陸海地形を含む
ようなケースでの大規模地震波シミュレーショ
ンを、東京大学FX10を用いてこれまでのプログラムよりも高速かつ安定に実行できることを確認できた。このことから、今後、陸海地形を含む構造モデルのもとでの沈み込み帯の巨大地震について、実用的な規模の計算を繰り返し行えることになる。

7．研究成果リスト

６．今年度の進捗状況と今後の展望(3) 国際会議発表
Taro Okamoto, Hiroshi Takenaka, Tatsuhiko Hara,
Takeshi Nakamura, and Takayuki Aoki, Rupture
Process of The 2011 Tohoku-Oki Earthquake
Inferred by an Inversion Using 3D Strong-Motion
and 2.5D Teleseismic Green's Tensor Waveforms,
AOGS 11th Annual Meeting, SE30-D1-AM1CA-001 (SE30-A016), Royton Sapporo Hotel,
Sapporo, Japan, 28 Jul. 2014.
Taro Okamoto, Hiroshi Takenaka, Takeshi Nakamura
and Takayuki Aoki, Large-Scale Multi-GPU
Simulation of Seismic-Wave Propagation for The
Shallow Subduction-Zone Earthquakes, AOGS 11th
Annual Meeting, SE36-D2-PM1-CA-003 (SE36A019), Royton Sapporo Hotel, Sapporo, Japan, 29
Jul. 2014.
Taro Okamoto, Hiroshi Takenaka, Tatsuhiko Hara,
Takeshi Nakamura, and Takayuki Aoki,
Seismically Inferred Rupture Process of The 2011
Tohoku-Oki Earthquake Inferred by Using DataValidated 3D and 2.5D Green's Tensor Waveforms,
American Geophysical Union 2014 Fall Meeting,
S33B-4514, San Francisco, California, USA, 17
Dec. 2014.

(1) 学術論文
Tanimoto, T. and T. Okamoto, The Millikan shaking
experiments and high-frequency seismic wave
propagation in Southern California, GeophysicalJournal International, 198, 1081-1095, doi:
10.1093/gji/ggu189.
(2) 国際会議プロシーディングス
該当なし。

• GPUカーネル部分の最適化（ 5.1）
各種の最適化について実測値に基づく検討を行い、最良のブロック個数・ブロックサイズ・
レジスタ数などに関する知見を得た。地震波計算プログラムのこれまでの実測最良性能は1GPUで226.6GFlopsである（K20X、CUDA6.0）。これは京コンピュータの1CPU理論ピー
ク性能のおよそ2倍であり、FX10の1CPU理論
ピーク性能にほぼ匹敵するものである。
• 短周期地震波再現性の検討（ 5.2）
小地震について観測波形と大規模シミュレー
ション結果とを用いた逆問題解析を行い、仮定した構造モデルのもとで観測波形を再現可
能な周期帯に関する知見を得た（周期約12秒
程度以上）。また、再現性向上のための構造
モデル改良の方向性についてある程度の見通
しを得ることができた。
• 摂動カーネル計算の効率化（ 5.3）
シミュレーションプログラムの出力において、
時刻歴波形ではなくて周波数領域に変換した
うえで少数の限定した周波数成分のみを出力することによりデータサイズを圧縮するようにプログラム改訂を行った。そのうえで大規模シミュレーションにより摂動カーネル計算を行い、結果について考察を行った。
• 地震波・津波統合計算（ 5.4）
地震波演算部を高速化するための最適化を施したプログラムを用いて、東京大学FX10で旧プログラムよりも高速化できることを確認した。この結果これまでは実行できなかった規模の計算がを容易に行うことができる見通しを得た。なお、津波項を含めた全体の高速化には至らなかったため、津波項を含む高速化の検討や津波項を含む大規模計算は次年度に継続して実施することとした。

(4) 国内会議発表
藤岡 慧・竹中博士・山田伸之・中村武史・岡元
太郎・藤原広行、南西諸島における強震動シ
ミュレーションのための３次元地下構造モデ
ル、 日 本 地 球 惑 星 科 学 連 合 2 0 1 4 年 大 会 、
SSS23-08、パシフィコ横浜（神奈川県横浜市）、2014年5月1日.
谷本俊郎・岡元太郎、アクティブソースによる高周波地震波の伝搬の研究、日本地球惑星科学
連合2014年大会、SSS27-14、パシフィコ横浜（神奈川県横浜市）、2014年4月29日.
岡元太郎・竹中博士・原辰彦・中村武史・青木
尊之、不均質構造モデルによる地震波グリーン関数を用いた2011年東北地方太平洋沖地震の破壊過程（4）、日本地震学会2014年秋季
大会、A31-01、朱雀メッセ：新潟コンベンショ
ンセンター（新潟市）、2014年11月26日.
岡元太郎・竹中博士・原辰彦・中村武史・青木
尊之、地震波から推定される2011年東北地方!8学際大規模情報基盤共同利用・共同研究拠点 平成26年度共同研究 最終報告書 2015年5月

太平洋沖地震の破壊過程：不均質構造モデルにもとづくグリーン関数を用いたインバージョ
ン解析、第14回日本地震工学シンポジウム、OS1-Sat-AM-5、幕張メッセ国際会議場（千葉市）、2014年12月6日.

(5) その他（特許，プレス発表，著書等）
該当なし。

参考文献
河井博紀、高速演算記 第25回「Kepler解説その2 〜Kepler世代の新機能〜」、日本GPUコン
ピューティングパートナーシップ（2013）［http://www.gdep.jp/column/view/33］

Okamoto, T. and H. Takenaka, Effect of near-source
trench structure on teleseismic body waveforms: anapplication of a 2.5D FDM to the Java trench, in
Advances in Geosciences, Vol.13 (Solid Earth), ed.
Kenji Satake, p.215–p.229, doi:
10.1142/9789812836182_0015, World ScientificPubishing Co., Singapore (2009).
Okamoto, T., H. Takenaka, T. Nakamura, and T.
Aoki, Accelerating large-scale simulation ofseismic wave propagation by multi-GPUs andthree-dimensional domain decomposition, Earth,Planets and Space, 62, 939-942 (2010).
Okamoto, T., H. Takenaka, T. Nakamura, and T.
Aoki, Accelerating large-scale simulation ofseismic wave propagation by multi-GPUs andthree-dimensional domain decomposition, in
GPU Solutions to Multi-scale Problems in
Science and Engineering, D.A. Yuen et al. (eds.),
Chapter 24, 375-389, doi:10.1007/978-3-64216405-7_24, Springer-Verlag Berlin Heidelberg(2013).
岡元太郎、大規模並列GPU計算による海溝型巨大地震の破壊過程推定と地震動評価、HPCI公募型共同研究成果報告書、2014.
Takenaka, H., Nakamura, T., Kuramoto, T., Okamoto,T., and Toyokuni, G., New FDTD Scheme for
Unified Modeling of Nearfield Tsunami and
Seismic Waves Using a Three-Dimensional EarthModel, APEC Cooperation for Earthquake
Simulation (ACES) 8 th International Workshop
Advances in Simulation of Multihazards, 9426,
Maui, Hawaii, Oct. 23, 2012.
Talantola, A., Inversion of seismic reflection data in
the acoustic approximation, Geophysics, 49, 12591266, 1984.
Tanimoto, T., Modelling curved surface wave paths:
membrane surface wave synthetics, Geophys. J.
Int., 102, 89-100, 1990.
Tromp, J., Tape C., and Liu, Q., Seismic tomography,
adjoint methods, time reversal and bananadoughnut kernels, Geophys. J. Int., 160, 195-216,2005.

!9