9:00 AM - 9:15 AM
[S15-10] Estimation of fault slip distribution using back propagation calculation of seismic waves by data assimilation - Numerical experiments
はじめに
一般に、大地震の震源断層すべり分布の推定は、地震波伝播のグリーン関数を用いたインバージョン解析に基づき行われる。これに対し、地震アレイ観測の波形データから震源断層すべり分布を直接イメージングするBack Projection（BP）が、近年、2004年スマトラ島地震（Ishii et al, 2005）や東北地方太平洋沖地震（Naercklin et al., 2012）などの巨大地震の震源過程解析に活用されるようになった。BPでは、仮定した断層面上の各点からの地震波放射に対応する見かけ速度でアレイ記録をスタッキングし、断層面上のエネルギー放射領域の時・空間変動を推定するものである。
さらに、アレイ観測記録を運動方程式の時間反転計算により震源に向けて発震時まで逆伝播させることで震源イメージングを行うTime reversal seismic source imaging （TRI）の近地・遠地アレイ観測データへの適用がLarmat et al. (2006)らによって提案されて以来、2004年スマトラ島地震（Larmat et al., 2006）や氷河の崩壊に伴う地震（Larmat et al., 2008）、火山性微動（Larmat et al., 2009）など各種の震源のイメージングへの活用を通して有効性が示されてきた。TRIでは、３次元地震波伝播シミュレーションに大きな計算資源を必要とするが、不均質な伝播経路の影響を考慮して、全波形記録の震源に向けた逆伝播を行うことで、震源イメージングの高解像度化が期待できる。しかし、Kremers et al. (2011)は各種数値実験に基づきTRIは点震源のイメージングには有効だが、その空間分解能、特に深さ方向の解像度に限界があり、大地震の地震断層すべり分布推定への適用は困難であると述べている。
データ同化を用いた新しいTRI
前述のTRIでは、地震波逆伝播計算にアレイ観測波形を直接または空間補間して入力される。このとき、観測点間隔が地震波の波長に比べて十分小さくないと球面状にばらばらに広がる波動場からは逆伝播波動場が適切に再構築されない。また、地震波逆伝播計算から求められる波動場は、観測点の位置で観測波形に完全に置き換えられるため、スムーズな波動場が再構築されにくい問題もある。
そこで本研究では、TRIにデータ同化の仕組みを導入し、逆伝播計算の結果と観測データが持つ誤差分布を考慮して適切に同化することでこの問題を解決し、逆伝播波動場を適切に再構築して震源イメージングの空間解像度を高めることを試みた。データ同化手法には、近年の津波（Maeda et al., 2015）や地震波（Hoshiba and Aoki, 2015; Furumura et al., 2019; Oba et al., 2020）の予測で広く用いられる最適内挿法を利用した。データ同化を用いたTRIは、既に古村・前田（2020）によって、中小地震の点震源のイメージングに適用され有効性が示されている。本研究では、これを大地震の断層すべり分布解析への適用に向け、その有効性を数値実験により検討した。
地震断層すべりイメージングの数値実験
まず、データ同化TRIの数値実験に用いる理論地震波形を２次元差分法計算により用意した。計算領域は、南海トラフを南北に横断する三重県沖〜日本海の650 km×50 kmの範囲とし、JIVSM（Koketsu et al., 2012）を用いて、堆積層・地殻・マントル構造とフィリピン海プレートを組み込んだモデルを用意した（図1）。設定した36 km×12 kmの震源断層モデルには２つのアスペリティを設定し、右下のアスペリティは、左上のおよそ9/10の大きさである（図１赤の長方形）。破壊開始点は断層中心に置き、破壊伝播速度Vr=3.2 km/sで断層面上の各格子点から地震波を放射させた。断層メカニズムは鉛直縦ずれ断層とし、地震発生から145秒間の地震波形を計算した。
次に、TRIの数値実験では時間を反転させた２次元差分法計算により、約10 km間隔の観測点から地震波を逆伝播させた。ここで、地震波動場の誤差分布が相関距離20 kmのガウス分布に従うと仮定し、また観測波形と計算結果の誤差の大きさが同程度であると仮定して、最適内挿法の重み関数を設定した。
データ同化を用いたTRI計算により得られた地震断層すべりイメージを図２に示す。ここでは、TRIの各時間ステップにおいて差分計算の各格子点の弾性エネルギーを求め、上位5格子点の合計エネルギーが最大となった時刻（4.5秒）での弾性エネルギーの分布を図示している。与えた断層すべりモデル（図１）と比較して、２つのアスペリティが分離した断層すべり分布のイメージングが適切に行われていることが確認できる。これに対して、データ同化を行わない、従来のTRIから求められた地震断層すべりイメージは、アスペリティを解像できずKremers et al. (2011)が指摘したとおりの結果となった。
なお、本実験では簡単のために、地震断層すべりの時間変化を考慮せず、弾性エネルギーが最大となった時刻の地震断層すべりを推定したが、巨大地震の震源解析では、断層すべりの時・空間的変化のイメージングが必要である。また、用いる地震波の周波数帯や観測点密度、またデータ同化パラメータがTRIの時間分解能に影響することが考えられる。数値実験によりこれらを確認し、大地震の実観測データを用いたTRIの有効性の検証を行う必要がある。
一般に、大地震の震源断層すべり分布の推定は、地震波伝播のグリーン関数を用いたインバージョン解析に基づき行われる。これに対し、地震アレイ観測の波形データから震源断層すべり分布を直接イメージングするBack Projection（BP）が、近年、2004年スマトラ島地震（Ishii et al, 2005）や東北地方太平洋沖地震（Naercklin et al., 2012）などの巨大地震の震源過程解析に活用されるようになった。BPでは、仮定した断層面上の各点からの地震波放射に対応する見かけ速度でアレイ記録をスタッキングし、断層面上のエネルギー放射領域の時・空間変動を推定するものである。
さらに、アレイ観測記録を運動方程式の時間反転計算により震源に向けて発震時まで逆伝播させることで震源イメージングを行うTime reversal seismic source imaging （TRI）の近地・遠地アレイ観測データへの適用がLarmat et al. (2006)らによって提案されて以来、2004年スマトラ島地震（Larmat et al., 2006）や氷河の崩壊に伴う地震（Larmat et al., 2008）、火山性微動（Larmat et al., 2009）など各種の震源のイメージングへの活用を通して有効性が示されてきた。TRIでは、３次元地震波伝播シミュレーションに大きな計算資源を必要とするが、不均質な伝播経路の影響を考慮して、全波形記録の震源に向けた逆伝播を行うことで、震源イメージングの高解像度化が期待できる。しかし、Kremers et al. (2011)は各種数値実験に基づきTRIは点震源のイメージングには有効だが、その空間分解能、特に深さ方向の解像度に限界があり、大地震の地震断層すべり分布推定への適用は困難であると述べている。
データ同化を用いた新しいTRI
前述のTRIでは、地震波逆伝播計算にアレイ観測波形を直接または空間補間して入力される。このとき、観測点間隔が地震波の波長に比べて十分小さくないと球面状にばらばらに広がる波動場からは逆伝播波動場が適切に再構築されない。また、地震波逆伝播計算から求められる波動場は、観測点の位置で観測波形に完全に置き換えられるため、スムーズな波動場が再構築されにくい問題もある。
そこで本研究では、TRIにデータ同化の仕組みを導入し、逆伝播計算の結果と観測データが持つ誤差分布を考慮して適切に同化することでこの問題を解決し、逆伝播波動場を適切に再構築して震源イメージングの空間解像度を高めることを試みた。データ同化手法には、近年の津波（Maeda et al., 2015）や地震波（Hoshiba and Aoki, 2015; Furumura et al., 2019; Oba et al., 2020）の予測で広く用いられる最適内挿法を利用した。データ同化を用いたTRIは、既に古村・前田（2020）によって、中小地震の点震源のイメージングに適用され有効性が示されている。本研究では、これを大地震の断層すべり分布解析への適用に向け、その有効性を数値実験により検討した。
地震断層すべりイメージングの数値実験
まず、データ同化TRIの数値実験に用いる理論地震波形を２次元差分法計算により用意した。計算領域は、南海トラフを南北に横断する三重県沖〜日本海の650 km×50 kmの範囲とし、JIVSM（Koketsu et al., 2012）を用いて、堆積層・地殻・マントル構造とフィリピン海プレートを組み込んだモデルを用意した（図1）。設定した36 km×12 kmの震源断層モデルには２つのアスペリティを設定し、右下のアスペリティは、左上のおよそ9/10の大きさである（図１赤の長方形）。破壊開始点は断層中心に置き、破壊伝播速度Vr=3.2 km/sで断層面上の各格子点から地震波を放射させた。断層メカニズムは鉛直縦ずれ断層とし、地震発生から145秒間の地震波形を計算した。
次に、TRIの数値実験では時間を反転させた２次元差分法計算により、約10 km間隔の観測点から地震波を逆伝播させた。ここで、地震波動場の誤差分布が相関距離20 kmのガウス分布に従うと仮定し、また観測波形と計算結果の誤差の大きさが同程度であると仮定して、最適内挿法の重み関数を設定した。
データ同化を用いたTRI計算により得られた地震断層すべりイメージを図２に示す。ここでは、TRIの各時間ステップにおいて差分計算の各格子点の弾性エネルギーを求め、上位5格子点の合計エネルギーが最大となった時刻（4.5秒）での弾性エネルギーの分布を図示している。与えた断層すべりモデル（図１）と比較して、２つのアスペリティが分離した断層すべり分布のイメージングが適切に行われていることが確認できる。これに対して、データ同化を行わない、従来のTRIから求められた地震断層すべりイメージは、アスペリティを解像できずKremers et al. (2011)が指摘したとおりの結果となった。
なお、本実験では簡単のために、地震断層すべりの時間変化を考慮せず、弾性エネルギーが最大となった時刻の地震断層すべりを推定したが、巨大地震の震源解析では、断層すべりの時・空間的変化のイメージングが必要である。また、用いる地震波の周波数帯や観測点密度、またデータ同化パラメータがTRIの時間分解能に影響することが考えられる。数値実験によりこれらを確認し、大地震の実観測データを用いたTRIの有効性の検証を行う必要がある。