第 7 章 「地殻活動シミュレーション手法」研究計画
１．はじめに
地殻活動シミュレーション研究の目標は，隣接するプレート同士が複雑に相互作用する日本列島域の地殻活動のシミュレーション・モデルを構築し，広域GPS観測網や地震観測網等からの膨大な地殻活動データをリアルタイムで解析・同化することにより，プレ−ト相対運動によって駆動されるテクトニック応力の蓄積から準静的な破壊核の形成を経て動的破壊の開始・伝播・停止に至る大地震発生過程の定量的な予測を行うことにある．
上記の目標を達成するためには，全国の大学及び関係諸機関が適切な役割分担の下に連携・協力し，複数の要素モデルをシステム結合した日本列島域を対象とするプロトタイプの地殻活動統合シミュレーション・モデルを構築する一方，大学等の研究グループが中心となって，以下に挙げるようなモデリング及びシミュレーション手法の高度化のための基礎研究を重点的に推進し，その成果を統合シミュレーション・モデルに逐次組み込むことで地殻活動予測シミュレーション・モデルを継続的に改良・発展させていく必要がある．
［大学等の研究グループが重点的に推進すべき基礎研究項目］
１）断層破砕帯の素過程
２）断層間相互作用
３）内陸活断層の地震発生過程
４）地殻活動データの解析・同化
５）特定地域に於ける地震発生サイクル・モデルの開発
６）日本列島域の広域変形・応力場のシミュレーション
２．平成13年度の研究成果
平成13年度の地震予知研究事業では，大学等の研究グループが重点的に推進すべきモデリング及びシミュレーション手法高度化のための基礎研究項目の内，１）の断層破砕帯の素過程に関する研究，２）の断層間相互作用に関する研究，３）の内陸活断層の地震発生過程に関する研究，5) の特定地域に於ける地震発生サイクル・モデルの開発，及び６）の日本列島域の広域変形・応力場のシミュレーションに関する研究が実施され，以下に示すような成果が得られた．
2.1 断層破砕帯の素過程に関する研究
平成12年度からの継続課題として「地殻内流体の挙動とその地震発生に対する力学的効果に関する研究（課題番号0131）」を実施した．平成13年度は，余震の発生に流体が及ぼす影響について数値シミュレーションを行い，余震発生の複雑さと規則性の生成要因に関して以下の成果を得た（Yamashita, 2002; 山下，2002）．
1)余震の発生に流体が及ぼす影響について数値シミュレーションを行い，大森公式およびグーテンベルグ・リヒターの式が統一的に再現できることを示した．
2)グーテンベルグ・リヒターの式を満たす地震は繰り返しすべりを起こしている破壊であることがわかった．
3)余震系列については，初期には比較的大きなイベントが起きる傾向があるなど，観測事実と調和的な結果が得られた．
4)複数の流体源がある場合，或いは未破壊領域の透水性がゼロに近いような場合は，二次余震が生じ得ることがわかった．
2.2 断層間相互作用に関する研究
平成12年度からの継続課題として「断層間相互作用による断層成熟度の変化についての研究（課題番号0132）」を実施する一方，平成11年度からの継続課題「地殻応力・歪変化シミュレーション手法に関する研究（課題番号0120）」の中で新たに断層間相互作用のシミュレーション研究を実施し，以下の成果を得た．
(1)断層間相互作用による断層成熟度の変化
同一平面上にはない複数の亀裂の動的な合体過程のシミュレーションを行い，一定速度で進む亀裂の相互作用による断層形状の変化を明らかにした．計算によると生成する断層形状は初期配置に依存し，二つの断層が離れ合う場合と近付き合う場合があることがわかった．また，断層形状は破壊伝播速度にも依存し，破壊伝播速度が大きくなると断層の屈曲が大きくなることもわかった（安藤・山下，2002）．
(2) 断層摩擦構成則に基づく地震サイクルシミュレーション
断層の速度・状態依存摩擦構成則を用いて，摩擦熱がプレスリップに及ぼす効果に関する研究（Kato, 2001a）や東海地震のプレスリップに関するシミュレーション（有吉・他, 2001）を実施した他，弾性層・粘弾性半空間中の横ずれ断層と無限媒質中の2次元断層面での地震サイクルの２種のモデルについて地震サイクルのシミュレーションを実施した．
1)弾性層とマックスウェル粘弾性半空間からなる2次元モデルを考え，鉛直横ずれ断層にはたらく摩擦力はすべり速度・状態依存摩擦則に従うと仮定して，地震サイクルシミュレーションを行った（Kato, 2001b）．断層の非地震性すべりと媒質の粘弾性的変形の両方を考慮することにより，インターサイスミック期の変形をより厳密にモデル化した．地震発生層の深さ下限が弾性層の厚さとほぼ同じである場合を除いて，断層でのすべりの時空間変化は弾性体半空間モデルでのものとほぼ同じであった．地震後の地表でのせん断歪変化に関しては，2つのモデルの違いは地震発生直後，粘弾性体の緩和時間に相当する期間にのみ顕著である．
2) 2次元断層面において，アスペリティのすべり過程や複数のアスペリティの相互作用を調べる目的で，均質無限媒質中のすべり速度・状態依存摩擦則に従う2次元平面断層での地震サイクルシミュレーションを行った（Kato et al., 2002）．定常すべり時の摩擦強度のすべり速度依存性を表すA-B，摩擦力のすべり量依存性を表すL，有効法線応力の空間的不均一を導入した．例としてA-Bの不均一を導入したシミュレーション（図1）の結果を示す．非地震性すべりをおこすA-B > 0の領域の中に，A-Bの絶対値が異なるA-B < 0の領域（アスペリティ）２つを埋め込んだ．A-Bの絶対値が小さいアスペリティはより早く非地震的にすべる．このすべりによる応力変化により，A-Bの絶対値が大きいアスペリティで地震が発生する．アスペリティで非地震性すべりがおこるか地震が起こるかは，アスペリティの大きさが摩擦パラメター値からきまる臨界サイズよりも大きいか小さいかで決まる（図2）．
2.3 内陸活断層の地震発生過程に関する研究
平成12年度からの継続課題として「下部地殻流動特性とプレート内応力の蓄積・解放過程のシミュレーション研究（課題番号0703）」を実施し，以下の成果を得た．
1)日本列島域の現実的な３次元プレート境界モデルを用いて，定常的なプレート沈み込みに伴う日本列島域の長期的地殻変動のシミュレーションを行い（図3），内陸活断層地震の原因となる地殻変形及び応力蓄積のパターンがプレート境界面の３次元形状に強く依存することを明らかにした（Hashimoto and Matsu’ura, 2002b）．また，このシミュレーションにより，既に行った横ずれプレート境界での地震発生サイクルのモデル化（Aochi and Matsu’ura, 2002; Hashimoto and Matsu’ura, 2002a; Fukuyama, Hashimoito and Matsu’ura, 2001）の場合と同様，プレート沈み込み帯の地震の発生サイクルを定常的なプレート沈み込み運動からの摂動分としてモデル化することが可能となった（Matsu’ura, 2001, 2002）．
2)内陸活断層の地震発生過程のモデル化に関しては，異方的な流動特性を持つ粘弾性物体の力学的応答の定式化とそれに基づく数値計算アルゴリズムの開発を進める一方，デタッチメントから分岐する衝上断層の地質学的時間スケールでの形状発達シミュレーションを行った（図4）．このシミュレーションにより，リストリックな内陸活断層の平均傾斜角が，デタッチメントでの定常すべりにより，低角から高角へと時間発展していくことが明らかになった（高田・松浦，2002）．
2.4 特定地域に於ける地震発生サイクル・モデルの開発に関する研究
平成13年度に課題内容を一部修正した「海溝型巨大地震の地震サイクルモデリング研究（課題番号0908）」を実施した．本年度は，海溝型巨大地震の準静的発生サイクルのモデルをめざして，東北日本における準静的な応力の蓄積過程に関する研究とブロックーバネモデルで単純化した南海トラフ沿いの巨大地震の発生シミュレーション研究を行い，以下の成果を得た．
(1) 東北日本における準静的な応力の蓄積過程
東北日本における応力蓄積レートの計算結果と観測された最近100年間の水平歪み速度との比較から，1900年頃三陸沖でモーメント・マグニチュード8.4に相当する大規模非地震性すべりが発生していたことが分かった（図5）．最近のGPS観測により，ゆっくり地震や地震後のアフタースリップの存在が報告されているが，このように大規模の非地震性のすべりは初めての発見である．またこれまで，三陸沖ではプレート間の固着は30％程度で，残りの70％は定常的にすべっていると考えられていた．しかし，観測された長期および短期の地殻変動を合理的に説明するには，プレート間の固着は100％で，その内の50〜60％が非地震性すべりを含むプレート境界地震で解放され，残りの40〜50％が大陸プレートを引きずり込むことで内陸地殻の応力蓄積に寄与しているとすればよい．内陸地殻の応力蓄積レートは，傾斜角45°の逆断層すべりを想定すると0.005〜0.01MPa/yr程度で，内陸大地震の応力降下量を10MPaとすると，その再来周期は数千年と見積もられ，地質学的に推定されている値とほぼ一致する（水藤尚・飯塚幹夫・平原和朗，2001）．
(2) ブロックーバネモデルによる南海トラフ沿いの巨大地震発生シミュレーション
南海トラフにそって発生する巨大地震の時空間的活動特性を理解するために，摩擦構成則に従った５つのセグメントから成るブロックーバネモデルの挙動を，パラメータマッチングの方法で調べた（光井能麻・平原和朗，2001, 2002）．図6aに示すようなブロックーバネモデルを設定し，同図右下の表にあるパラメータを用いて地震サイクルのシミュレーションを行った．その結果，ブロックA〜Dがほぼ同時にすべる場合（図6b左）が2/3の頻度で発生し，同時ではなく交互にすべる場合（図6b右）が１/３の頻度で発生することが分かった．また，前者の場合はその繰り返し間隔が60年〜120年程度のばらつきを持つこと（図6c左），後者の場合は BとCの地震発生時間差の80％が1年以内であること（図6c右）が示された．南海トラフでのデータは1,000年程度なので，シミュレーション結果のようなことが実際生じるのかどうか定かではないが，大きく異なる地震発生サイクル状態が存在する可能性がある．
2.5 日本列島域の広域変形・応力場のシミュレーションに関する研究
平成11年度からの継続課題として「地殻応力・歪変化シミュレーション手法に関する研究（課題番号0120）」を実施し，平成13年度は逆解析手法に基づく地殻応力変化の推定に関して以下の成果を得た．
ＧＰＳ観測から得られる変位場から，逆解析手法を用いて，応力変化をモニタリングするシステムを構築した（Hori et al., 2001; Iinuma, 2001）．この手法では地殻を薄い弾塑性媒質で近似し，平面応力場を仮定する．変形に関するつりあいの条件と塑性変形部に関する構成関係を仮定することにより問題を定式化し，Airyの応力関数を導入することでGPSデータから精度良く地殻応力を推定することが可能となった．1996〜1999の３年間のGPSデータを用い，日本列島の応力増分を推定した．推定された応力がどの程度確からしいかを判定するため，推定された応力から歪みを逆算し，弾性体を仮定してGPSデータから直接計算した歪みとの差をとった（こうすることで地殻の剛性率の分布を間接的に推定できる）．これを同じ期間に発生した地殻内地震と比較した（図7）．その結果，剛性率が小さい場所で地震活動が高い傾向にあることが見て取れた．
３．まとめ
平成12年度から「地殻活動シミュレーション手法」に関する三つの新たな研究課題（0127, 0128, 0703）が地震予知研究事業として追加されたことにより，平成13年度は，モデリングやシミュレーション手法の高度化のための基礎研究が活性化された．一方，全国の大学及び関係諸機関の研究者の連携・協力の下に，科学技術振興調整費総合研究（平成10〜15年度）として開発中の日本列島域を対象とする「地殻活動統合シミュレーション・モデル」のプロトタイプは，平成13年度末の地球シミュレータの完成により，平成14年度には組み上がる予定である．上記の平成13年度の研究成果は，平成14年度以降の継続・進展分も含めて，このプロトタイプの統合シミュレーション・モデルを改良・発展させていくための重要な基礎研究として位置付けることができる．
参考文献：
Aochi, H. and M. Matsu'ura, Slip- and time-dependent fault constitutive law and its significance in earthquake generation cycles, Pure Appl. Geophys. 159, 2029-2047, 2002.
有吉慶介・加藤尚之・長谷川昭，東海地域における近年の地殻変動及び地震活動の変化に関する数値シミュレーションによる検討，地学雑誌，110, 557-565, 2001．
安藤亮輔・山下輝夫，断層の幾何学的形状の形成と動的破壊過程，地球惑星科学関連学会2002合同大会，東京（2002年5月27日—31日），S040-010, 2002．
Fukuyama, E., C. Hashimoto and M. Matsu'ura, Simulation of earthquake rupture transition from quasi-static growth to dynamic propagation, Proceedings of the 2nd ACES Workshop, ed. M. Matsu’ura, K. Nakajima and P. Mora, APEC Cooperation for Earthquake Simulation, 375-380, 2001.
Hasimoto, C., and M. Matsu'ura, 3-D simulation of earthquake generation cycles and evolution of fault consitutive properties, Pure Appl. Geophys. 159, 2175-2199, 2002a.
Hasimoto, C., and M. Matsu'ura, Long-term crustal deformation in and around Japan, simulated by a 3-D plate subduction model, Abstracts of the 3rd ACES International Workshop, Maui (May 5-10, 2002), APEC Cooperation for Earthquake Simulation, 2002b.
Hori, M., T. Kameda, and T. Kato, Application of the inversion method to a GPS network for estimating the stress increment in Japan, Geophys. J. Int., 144, 597-608, 2001.
Iinuma, T., Stress inversion analysis based on the velocity field of the Japanese Islands, Ms. Thesis, 2002.
Kato, N., Effect of frictional heating on preseismic sliding: A numerical simulation using a rate, state, and temperature dependent friction law, Geophys. J. Int., 147, 183-188, 2001a.
Kato, N., Seismic cycle on a strike-slip fault with rate- and state-dependent strength in an elastic layer overlying a viscoelastic half-space, Proc. Int. Symp. Slip and Flow Processes in and below the Seismogenic Region, 103-108, 2001b.
Kato, N., K. Hirahara, and M. Iizuka, Numerical simulation of seismic cycles at a subduction zone with a rate- and state-dependent friction law, Abstracts of the 3rd ACES International Workshop, Maui (May 5-10, 2002), 37, 2002.
Matsu'ura, M., The crustal activity modelling program: Progress toward scientific forecast of earthquake generation, Proceedings of the 2nd ACES Workshop, ed. M. Matsu’ura, K. Nakajima and P. Mora, APEC Cooperation for Earthquake Simulation, 23-26, 2001.
Matsu’ura, M., Development of physics-based predictive simulation model for earthquake generation, Abstracts of the 3rd ACES International Workshop, Maui (May 5-10, 2002), 41, 2002.
光井能麻・平原和朗，沈み込み帯における地震発生サイクルシミュレーションーコントロールオパラメータとそれに対応したモデルパラメータの解釈—，日本地震学界講演予稿集，2001年度秋季大会，A46, 2001．
光井能麻・平原和朗，バネ−ブロックモデルによる南海トラフ巨大地震発生サイクルシミュレーション，号外地球，2002（印刷中）．
水藤尚・飯塚幹夫・平原和朗，東北日本の地殻変動シミュレーション，日本地震学界講演予稿集，2001年度秋季大会，B71, 2001．
高田陽一郎・松浦充宏，プレート境界面から分岐する衝上断層の時間発展：ヒマラヤのテクトニクスへの寄与，地球惑星科学関連学会2002年合同大会，東京（2002年5月27日—31日），T043-021, 2002．
Yamashita, T., Mechanical effects of fluid migration in a fault zone on seismic activity, Abstracts of International Workshop on Physics of Active Fault, Tsukuba (26-27 Feb., 2002), S3-5, 2002.
山下輝夫，断層帯内流体移動の地震活動に対する力学的効果，地球惑星科学関連学会2002合同大会，東京（2002年5月27日—31日），G061-018, 2002．
図の説明：
図１．2次元断層面上の摩擦パラメターA-Bの分布．
図２．1サイクルでの地震性すべり量の空間分布．
図３．定常的なプレート沈み込み運動が引き起こす日本列島域の長期的地殻変動のシミュレーション．(a) 日本列島域の３次元プレート形状モデル．(b) プレート沈み込みモデルで計算される長期的な地殻の隆起（赤）・沈降（青）速度．
図４．プレート収束運動によって駆動される分岐衝上断層の形状発達シミュレーション．(a) 分岐衝上断層を持つデタッチメントの定常すべりによる内部変位速度場（上磐側ブロックの無限遠点を固定）．(b) プレート収束に伴う分岐衝上断層の形状変化．
図５．仮想的な地震を含む100年間の歪速度の主軸の計算値（左）と三陸沖のモーメント解放の時空間分布（右）．
図 ６．ブロックーバネモデルを用いた南海トラフ巨大地震発生のシミュレーション． (a)モデルの設定（左；南海トラフ巨大地震発生サイクルとセグメント，右上；ブロックーバネモデル，右下；各セグメントのパラメータ）． (b)各ブロックにおける地震発生の時系列． (c)地震発生間隔と発生時間差の分布（左；５ブロックモデルの繰り返し時間間隔の度数分布，右；BとCの地震発生時間差の度数分布）．
図７．応力逆解析から推定した地殻ひずみと弾性体近似によるひずみとの差分（剛性率の不均質性）と地震活動（D<=20km; M>=2.0）．1996.04-1999.08のデータを使用．
以上
