第 3 章 手法
3.2 震源モデル解析
48
津波波源インヴァージョンおよび一枚矩形断層グリッドサーチの際には,海底の上下変 動分布をそのまま初期の海面の波高分布として津波計算を行う.また,断層すべりインヴ ァージョンの際には,海底面の上下変動場に対して,水深フィルター (Saito and Furumura,
2009) を適用し,それを初期海面波高分布として津波計算を行った. 詳細は次節において
述べる.
49
ンにより推定される初期波高分布が,既知の発震機構解から期待される分布になるか比較 することにより,分岐断層の破壊や海底地すべりなど,断層すべりでは説明できないよう な海底変動が起こっていないか検証できる,そこで,予備解析の一種として津波波源イン ヴァージョンを行い,その結果を踏まえて震源モデルの推定を行う.
インヴァージョンには,津波波形がグリーン関数の線形結合によって表現されることを 利用する.図3.4に示すように,観測される津波波形は,波形計算の単位要素,すなわち津 波波源インヴァージョンの場合には解析領域を小さく区分けした小波源が変動することに よって発生する津波の重ね合わせで表現される.このことは次式で表される観測方程式で 表現できる.
𝑑𝑗(𝑡) =� 𝐺𝑖𝑗(𝑡)𝑎𝑖
M
𝑖=1
(3.9)
ここで,
𝑑𝑗(𝑡):j番目の観測点における津波波形時系列データ
𝐺𝑖𝑗(𝑡):i番目の小波源による,j番目の観測点における津波波形時系列(グリーン 関数)
𝑎𝑖:i番目の小波源における変動量 M:小波源の総数
である.これを行列表示で書くと,次式で表される.
𝐝=𝐆𝐚 (3.10)
ここで,
N =�時間ステップ数�×�観測点数�
とすると
50 𝒅:𝑑𝑗(𝑡)を成分とするN行1列の列ベクトル 𝐆:𝐺𝑖𝑗(𝑡)を要素とするN行M列の行列 𝐚:𝑎𝑖を要素とするM行1列の列ベクトル
である.インヴァージョンを安定化させるために,空間的な平滑化の拘束を与える.具体 的には,隣り合う海底面要素の上下変位量のラプラシアンがゼロになるような拘束であり,
次式の形で表現される.
0 =𝑎𝑙+1,𝑚+𝑎𝑙−1,𝑚+𝑎𝑙,𝑚+1+𝑎𝑙,𝑚−1−4𝑎𝑙,𝑚 (3.11)
ここで,
𝑎𝑙,𝑚:x方向にl番目.y方向にm番目の海底面要素における上下変動量
である.行列表示で書くと,次式で表される.
𝟎=𝐒𝐚 (3.12)
ここで,Sは(3.11)式の各係数を要素とするN行N列の行列である.
平滑化の拘束条件を考慮した津波波源インヴァージョンの目的関数𝑠(𝐚)は,次式で表され る.
𝑠(𝐚) =‖𝐝 − 𝐆𝐚‖2+𝛼2‖𝐒𝐚‖2 (3.13)
ここで,𝛼は平滑化の拘束条件の重みを表す.この目的関数𝑠(𝐚)が最小の値をとるような最 適解𝐚∗を推定する.正規方程式は次式で表される.
(𝐆T𝐆+𝛼2𝐒T𝐒)𝐚=𝐆T𝐝 (3.14)
この式から,解a*を,特異値分解法により求める.𝛼の値は,本研究で解析を行う3つの地 震すべてについて経験的に求めた 𝛼= 1 を用いた.具体的には,𝛼の値を 0.1,1,10,… と変化させてインヴァージョンを行い,3つの地震すべてで解析領域の端にインヴァージョ
51
ン上の虚像が現れない値として,経験的に得た 𝛼= 1 を採用した.
続いて,津波波源インヴァージョンに用いたグリーン関数の計算の手法について述べる.
グリーン関数計算の元となる,単位源としての小波源には,図3.5に示した,隣り合う小断 層どうしがその大きさの半分ずつ重なるピラミッドのような形を用いた.すなわち,海底 面要素を20km四方の正方形とし,単位変位量には中心を1m,中心から離れるにつれて,
次式で表される変位量を与えた.
𝜂0(𝑥,𝑦) =�1− �𝑥
10�� �1− �𝑦 10��
ただし−10 <𝑥< 10,−10 <𝑦< 10
(3.15)
実際の津波の波源としての海底上下変動場は,上下変動量が水平方向になめらかに変化 するものと考えて,単純な直方体ではなくこのような形とした.この小波源を初期条件と して津波の数値計算を行う.
地震が発生した際,震源域の直上に設置された海底圧力計では,断層運動に伴う海底上 下変位が,津波による海面高変動と共に記録されることになり,図3.6のように地震前後で 水圧オフセットが変化する.本研究におけるグリーン関数の計算では,Tsushima et al. (2012) の手法に基づき,地震時の海底の上下変動に対応する圧力変化を,計算された津波から差 し引くようにしている.これにより,水圧オフセットの変化分を別途推定することなく,
津波波源インヴァージョンを行うことができる.
グリーン関数を計算する領域は地震によって異なるが,震源およびその周囲の観測点を 含む十分広い領域を対象とすることとした.
3.2.2 断層すべりインヴァージョン
本節では,断層のすべり分布を推定する断層すべりインヴァージョンについて述べる.
断層すべりインヴァージョンも初期波源インヴァージョンと同様に,津波波形がグリー
52
ン関数の線形結合によって表現されることを利用する.図3.4に示すように,観測される津 波波形は,波形計算の単位要素,断層すべりインヴァージョンの場合には断層面を小さく 区分けした小断層が変動することによって発生する津波の重ね合わせで表現される.断層 すべりインヴァージョンにおける観測方程式は,(3.9)式と同様の形で表される.
𝑑𝑘(𝑡) =� 𝐺𝑖𝑘(𝑡)𝑎𝑖
M
𝑖=1
(3.16)
ここで,
𝑑𝑘(𝑡):k番目の観測点における津波波形データ
𝐺𝑖𝑘(𝑡):i番目の小断層による,k番目の観測点における津波波形時系列 (グリーン 関数)
𝑎𝑖:i番目の小断層におけるすべり量 M:小断層の総数
である.これを行列表示で書くと,次式で表される.
𝐝=𝐆𝐚 (3.17)
ここで,
N =�時間ステップ数�×�観測点数�
とすると
𝐝:𝑑𝑘(𝑡)を成分とするN行1列の列ベクトル 𝐆:𝐺𝑖𝑗(𝑡)を要素とするN行M列の行列 𝐚:𝑎𝑖を要素とするM行1列の列ベクトル
である. 断層すべりインヴァージョンにおいても,津波波源インヴァージョンと同様に空 間的な平滑化の拘束を与える.その与え方は(3.11)式と同様である. 以上より,正規方程式
は,(3.14)式と同様に次式で表される.
53
(𝐆T𝐆+𝛼2𝐒T𝐒)𝐚=𝐆T𝐝 (3.18)
この式から,解 a*を,特異値分解法により求める.なお,断層すべりインヴァージョンで は,スムージングの重み𝛼は経験的に決めるのではなく,𝛼および計算波形と観測波形の残 差のトレードオフ曲線を描くことにより決定した.詳細は第4章で述べる.
断層のすべり分布やその広がりについて議論を行うために,断層すべりインヴァージョ ンにより求まる断層すべり量の推定誤差を評価する必要がある.i番目の小断層におけるす べり量𝑎𝑖の推定誤差σ𝑖𝐚は,観測データベクトル𝐝が,真の値からなるベクトル𝐝0と,期待値 0,分散共分散行列Σの正規分布に従う観測誤差ベクトル𝛜を用いて,
𝐝=𝐝𝟎+𝛜 (3.19)
と表されるとき,
𝜎𝑖𝐚= {[(𝐆T𝚺−𝟏𝐆+𝛼2𝐒T𝐒)−1]𝑖𝑖}12 (3.20)
と求められる.本研究では分散共分散行列Σの非対角成分を 0,すなわち,すべての観測誤 差が独立であるものとして,
𝚺= Cov[𝛜𝛜T] =
⎣⎢
⎢⎡𝜎12 𝐎 𝜎22
𝐎 ⋱ 𝜎N2⎦⎥⎥⎤
(3.21)
とした.ここで,𝜎𝑚2はm番目の観測データの持つ誤差の分散である.さらに,すべての観 測点,すべての時刻において観測誤差の分散が等しい,すなわち,
𝜎1=𝜎2=⋯=𝜎N=𝜎0 (3.22)
と仮定し,𝜎0= 1 hPa≅1 cmとした.
続いて,断層すべりインヴァージョンにおける,グリーン関数の計算について述べる.
断層のすべりインヴァージョンのグリーン関数の計算は,図3.4で表すように,津波波源イ ンヴァージョンの場合と異なり,断層運動による津波を考える必要がある.そこで,小断
54
層ごとに,そこでの単位すべりがもたらす海底の上下変動をOkada (1992) の式により計算 し,得られた初期海底上下変動場から津波を計算する.なお,本研究において断層すべり インヴァージョンを適用するのはイベント0309 およびイベント 0310 である.これらの地 震が発生した海域の水深は平均して2 kmほどであるので,Okada (1992) の式を用いて海底 の上下変動を計算する際に,自由表面は海面から2 kmの深さにあるとした.
Okada (1992) の式から海底上下変動場𝑢𝑧(𝑥,𝑦)を計算したのち,水深フィルター (Saito and
Furumura, 2009) を適用し,海面の初期波高分布を計算する.一般に,津波計算では初期波
高分布は,地震時の海底上下変動場と等しいと仮定することが多い.しかし,水深に対し て海底上下変動域の広がりが狭い場合はこの仮定が成り立たない.Saito and Furumura (2009) では,隆起域,または沈降域の広がりが,水深の13倍よりも狭いときには水深フィルター による初期波分高布の補正が必要になると述べている.水深フィルターは,次式で表され る.
𝜂(𝑥,𝑦) = F−1� 1
cosh𝑘ℎ 𝑢��𝑘𝑥,𝑘𝑦�� (3.23) ここで,𝑢��𝑘𝑥,𝑘𝑦�は𝑢𝑧(𝑥,𝑦)の空間でのフーリエ変換,すなわち,𝑢��𝑘𝑥,𝑘𝑦�= F[𝑢𝑧(𝑥,𝑦)]で あり,
𝑢��𝑘𝑥,𝑘𝑦�=� �∞d𝑥d𝑦 𝑢𝑧
−∞
∞
−∞ (𝑥,𝑦)exp�−𝑖�𝑘𝑥𝑥+𝑘𝑦𝑦�� (3.24) で表される.また,
𝒌=�𝑘𝑥,𝑘𝑦�: 波数ベクトル ℎ:水深 (一定)
である.本研究では,海底の上下変動の計算の時と同様に ℎ= 2 kmとした.このようにし て得られた初期海面変動場を用いて津波計算を行い,計算により得られた各観測点の時系 列から,海底の上下変動量に伴う圧力変化を差し引いたものをグリーン関数とした.
55
3.2.3 一枚矩形断層グリッドサーチ
本節では,一枚矩形断層グリッドサーチについて述べる.イベント0309やイベント0310 のように,プレート境界を震源として発生した逆断層型地震による津波波形を解析する場 合には,多くの先行研究(例えばIto et al., 2004;Yamamoto et al., 2011)によってプレート 境界面の形状があきらかとなっていることから,断層面形状としてプレート境界面の形状 を仮定し,その断層面上でのすべり分布を推定することが可能である.しかし,イベント 0710 は沈み込む太平洋スラブ内で発生した横ずれ地震であり,その震源となった断層のパ ラメタを十分な精度で拘束できるような先験情報が少ない.そこで,この地震については 一枚の平面矩形断層を用いたグリッドサーチにより震源モデルの推定を行い,断層の面積 だけでなく傾斜やすべり角も推定することとした.推定する一枚矩形断層の模式図および パラメタを,図 3.7 に示す.なお,3.2.1 で述べた津波波源インヴァージョンの結果に基づ いて,断層の走向の値は固定し,Global CMT の値を用いた.詳細は第5章にて述べる.
本研究では,計算の時間を節約するために,グリッドサーチを 2 つの段階に分けて行う こととした.まず 1 つ目の段階で,探索の範囲を広く,かつ粗く計算して探索の範囲を絞 り込み,2つ目の段階において,絞り込まれた範囲の中で詳細に探索を行う.
1 段階目のグリッドサーチでは,図 3.7 で表されるパラメタのうち,断層の走向は
GlobalCMT の値に固定し,断層の傾斜,すべり角,断層の長さ及び幅を変化させて矩形断
層を仮定し,Okada (1992) のモデルから海底上下変動を計算する.その際,自由表面がイ ベント0710の震源域の平均水深である3.5 kmの深さにあるものとして,計算を行ってい る.その後,次式で表される,計算した海底上下変動分布と津波波源インヴァージョンに よって得られた初期津波波高分布との残差二乗平均平方根 (Rooted Mean Square;RMS) を 計算する.