ＩＳＳＮ１３４６－７８４０

港湾空港技術研究所
資料
TECHNICAL NOTEOFTHE PORT AND AIRPORT RESEARCH INSTITUTE

No.1401

March 2022

フリングステップを考慮できる複合型震源モデルの提案
野津 厚，呉国立研究開発法人

双蘭，長坂

陽介

海上・港湾・航空技術研究所

National Institute of Maritime,
Port and Aviation Technology, Japan目

次旨 ························································································································

3

1. まえがき ···················································································································

4

2. 複合型震源モデル ·········································································································

52.1

複合型震源モデルの全体像 ······················································································

5

2.2

離散化波数法の概要 ·······························································································

6

3. 離散化波数法の計算条件 ·······························································································

7

要

3.1

永久変位に関する静的な解析解との比較 ····································································· 7

3.2

増田・引間による解析解との比較 ············································································ 11

4. 熊本地震本震への適用 ·································································································· 13
4.1

国土地理院による震源断層モデルとその適用 ····························································· 13

4.2

アスペリティの考慮 ······························································································ 17

5. まとめ ······················································································································· 17
6. おわりに ··················································································································· 21
参考文献 ························································································································ 21-1-

A Composite Source Model Capable of Generating Fling Steps
Atsushi NOZU*
Shuanglan WU**
Yosuke NAGASAKA***

Synopsis
During large crustal earthquakes, strong ground motions with a significant amount of residual
displacement are observed if the rupture of the fault approaches very close to or reaches the earth’s
surface. Those ground motions are often called ‘fling steps’. Fling steps are mainly due to the ruptureof shallower parts of the fault. During the mainshock of the 2016 Kumamoto earthquake, significant
fling steps were observed because the rupture reached the surface. On the other hand, during the 1995
Kobe earthquake, the rupture did not reach the surface in the Kobe side and the devastating groundmotions were mainly associated with the rupture of deeper asperities. This was the reason why theconventional studies related to design ground motions were mainly focused on the rupture of deeper
asperities and less attention was paid to fling steps. In this study, we propose a composite source
model which can consider fling steps in addition to strong ground motions from deeper asperities.
The source model comprises of Domain A representing deeper asperities and Domain B whichreaches the surface and generates fling steps. The ground motions from Domain A are calculated by
the corrected Empirical Green’s function method, while the ground motions from Domain B are
calculated by the discrete wavenumber method because the reflection coefficients at the free surface
should rigorously be considered for the calculation of the fling steps. In this study, first, the resultsof the discrete wavenumber method were compared with analytical results including Okada’s
analytical solutions for static crustal deformations to understand the criteria for accurately calculating
fling steps with the discrete wavenumber method. Then a composite source model was actuallydeveloped for the mainshock of the 2016 Kumamoto earthquake. The authors’ intention was todevelop a simple model which can be easily handled by engineers. For this reason, the geodetic fault
model by the GSI was used as a starting point and parameters related to the temporal evolution of theslip were added to form Domain B. Then it was combined with Domain A representing asperities toform a composite source model. The simulation results showed that the composite source model canwell reproduce near-source displacement and velocity waveforms including fling steps, although themodel was a relatively simple one in a sense that the slip-velocity function and the rise time were
uniform for each of the faults in Domain B.
Key Words: Fling step, source model, fault, earthquake ground motion, discrete wavenumber
method* Director of Earthquake Disaster Prevention Engineering Department** Postdoctoral Researcher, Engineering Seismology Group, Earthquake Disaster Prevention Engineering Department*** Senior Researcher, Engineering Seismology Group, Earthquake Disaster Prevention Engineering Department
3-1-1 Nagase, Yokosuka, 239-0826 Japan
Phone：+81-46-844-5058
Fax：+81-46-844-0839
e-mail:nozu@p.mpat.go.jp-2-フリングステップを考慮できる複合型震源モデルの提案
野津

厚*・呉

双蘭**・長坂

要

陽介***

旨

規模の大きい内陸地殻内地震において，断層の破壊が地表付近まで到達する場合や，地表地震断層が現れる場合は，その近傍において永久変位成分を含む地震動が生じる．このような永久変位成分を含む地震動はフリングステップと呼ばれることが多い．フリングステップは主に断層浅部のすべりにより生じるものであり，断層の破壊が地表に到達した 2016 年熊本地震本震では顕著であったが，断層の破壊が地表に到達しなかった 1995 年兵庫県南部地震の神戸側では，断層深部に位置するアスペ
リティの破壊による影響が支配的であり，フリングステップは顕著でなかった．このため，従来の設計地震動に関する研究では，深部に位置するアスペリティの破壊による影響が主に考慮され，浅部のすべりの影響は十分に考慮されてこなかった．そこで本研究では，従来から考慮されているアスペリ
ティの破壊に伴う強震動に加え，フリングステップも考慮できる複合型の震源モデルを新たに考案した．この震源モデルは，アスペリティに対応する領域 A とその外側の領域 B からなる．領域 A は深部に留まるが，領域 B は地表付近まで到達しフリングステップを生じる．領域 A からの地震動は従来通り修正経験的グリーン関数法で計算するが，領域 B からの地震動は離散化波数法で計算し，その結果を加算する．領域 B からの地震動を離散化波数法で計算するのは，フリングステップの計算において，地表面における反射係数を厳密に扱うことが必要だからである．本研究では，まず，離散化波数法による計算結果を静的な永久変位に関する解析解等と比較することで，フリングステップを正しく計算するために必要な計算条件を明らかにした．次に，複合型の震源モデルを実地震である熊本地
震本震に適用した．震源モデルをできるだけ簡便で利用しやすいものとすることを念頭におき，測地学的データを説明するように設定された国土地理院の震源モデルを出発点とし，これに最小限の動的な情報を付与することで領域 B を設定し，これに領域 A を加えることで複合型の震源モデルを構築したところ，実際の観測事実と整合するシミュレーション結果を得ることができた．
キーワード：フリングステップ，震源モデル，断層，地震動，離散化波数法* 地震防災研究領域長** 地震防災研究領域 地震動研究グループ 専任研究員*** 地震防災研究領域 地震動研究グループ 主任研究官〒239-0826 横須賀市長瀬3-1-1 国立研究開発法人海上・港湾・航空技術研究所 港湾空港技術研究所
電話：046-844-5058 Fax：046-844-0839
e-mail:nozu@p.mpat.go.jp-3-

1. まえがき

上述のように，社会基盤施設に大被害をもたらした1995年兵庫県南部地震の神戸市内の地震動に対しては，規模の大きい内陸地殻内地震において，断層の破壊が断層深部のアスペリティの破壊による影響が支配的であ

地表付近まで到達する場合や，地表地震断層が現れる場

ったことから，これまでの設計地震動に関する研究は，断合は，その近傍において永久変位成分を含む地震動が生層深部のアスペリティの破壊による影響を主に考慮する

じる．このような永久変位成分を含む地震動はフリングという形で進められてきており，断層浅部のすべりの影

ステップと呼ばれることが多い（Hisada and Bielak, 2003）．

響は十分に考慮されてこなかった面がある．そのため，現フリングステップは主に断層浅部のすべりにより生じ時点で一般に用いられている強震動シミュレーション手るものである．1995年兵庫県南部地震の際には，神戸側の法は，フリングステップを考慮できないものが多い．例え

断層は淡路側の断層ほど浅部のすべりは大きくなかった

ば，港湾の分野で広く照査用地震動の評価に用いられてため（例えばSekiguchi et al., 2000），神戸市内で観測されいる修正経験的グリーン関数法（古和田他，1998；野津・た地震動には顕著なフリングステップは含まれておらず，

菅野，2008；野津他，2009）は，兵庫県南部地震の神戸市

神戸市内の地震動に対しては，深部（数km程度以深）の内の地震動など，永久変位成分が支配的でない場合の震アスペリティの破壊による影響が支配的であった（例え

源近傍地震動の計算おいて実績があるが，そのままでは

ば入倉，1996）．しかしながら，2016年熊本地震本震の際

フリングステップを含む地震動を計算することができなには，地表まで破壊が進展したことに伴い（例えばい．こうした状況は，フリングステップの構造物への影響

Shirahama et al., 2016），断層近傍の西原村小森などでは顕を調べる必要性が増していることを考慮すると不都合で著なフリングステップが観測された（図-1）．

あると考えられる．

フリングステップは，その立ち上がり時間が十分に大そこで本研究では，従来から設計において考慮されてきければ，断層近傍に立地する構造物への影響は（断層をいるアスペリティの破壊に伴う強震動に加え，フリング

またぐ構造物の場合を除けば）小さいと考えられる．しかステップも考慮できる複合型の震源モデルを新たに考案しながら，熊本地震本震の際に観測されたフリングステ

し，これにより強震動シミュレーションを行うことを提ップは立ち上がり時間が2s程度と短く，断層をまたがない

案する．提案する震源モデルは，アスペリティに対応する

構造物に対しても顕著な応答をもたらしうるものであっ

領域Aとその外側の領域Bからなる．領域Aは深部に留また．こうしたことから，フリングステップの重要性が再認るが，領域Bは地表付近まで到達しフリングステップを生識されており，特に建築分野において，免震構造物等へのじる．領域Aからの地震動は従来通り修正経験的グリーン関数法で計算するが，領域Bからの地震動は離散化波数法（Bouchon, 1981; Bouchon, 2003）で計算し，その結果を加に作用させる研究が行われるようになってきている（例算する．領域Bからの地震動を離散化波数法で計算するのえば豊岡他，2021）．

は，フリングステップの計算において，地表面における反Displacement (cm)

影響が精力的に研究されている（例えば林・杉野，2018）．
土木分野においても，フリングステップを構造物に動的

150
100

Komori, Nishihara, Kumamoto Pref.

50

Integrated in time domain
Integrated in frequency domain

0
-50
0

10

20

30

40

50

60

Time (s)図-1 熊本地震本震の際に西原村小森で得られた強震記録を積分して得られた変位波形（野津，2018a），縦軸：変位，横軸：時間，灰色線：時間領域での積分結果，赤線：周波数領域での積分結果．この場合150cm程度の残留変位が生じていることがわかる．すなわち，揺れが収まった時には，地震前とは異なる位置で地盤が停止している．このような永久変位成分を含む地震動をフリングステップと呼ぶ．

-4-

本的考え方を示す．この震源モデルは，アスペリティに対射係数を厳密に扱うことが必要だからである．

応する領域Aとその外側の領域Bからなる．領域Aは深部

以下においては，まず，離散化波数法による計算結果を
静的な永久変位に関する解析解等と比較することで，フに留まるが，領域Bは地表付近まで到達しフリングステッ

リングステップを正しく計算するために必要な計算条件

プを生じる．ここで言う深部とは少なくとも数km以深とを明らかにする．次に，複合型の震源モデルを実地震であいう意味である．もともと内陸地殻内地震は深さ20km程る熊本地震本震に適用し，その妥当性を検討する．

度の範囲で起きるが，耐震設計用の地震動を計算する際

一般に強震動シミュレーションの研究では，実現象のは，従来は地表面から3km程度の範囲は地震波を出さない再現を目指すのか，将来の予測や設計のために簡便で利と考え，かつ，領域Aからの地震動を計算するだけで足り用しやすい手法の開発を目指すのかによって，アプロー

ると考えられていた．実際，兵庫県南部地震などの揺れはチがかなり変わってくる．熊本地震本震のフリングステ

領域Aだけで十分に説明できてしまう．しかし，領域Aかップについても，実現象としてどのようなものであった

らの地震動は残留変位を説明できず，残留変位成分を顕かについての研究は他の研究機関によりかなり行われて著に含む地震動を説明するためには領域Bも必要となる．

いる．それらの研究では，断層面の形状もより複雑なもの

領域Bは深部も含め設定されるが，特に断層近傍での地震を考え，断層面上のすべり量分布についてもより複雑な動を計算する場合，浅部の影響が支配的となる．

ものを考えることで，実現象の再現を目指すといったア領域Aからの地震動は従来通り修正経験的グリーン関

プローチがとられる．しかし，仮に複雑な断層形状やすべ数法（古和田他，1998；野津・菅野，2008；野津他，2009）

り分布を有するモデルにより観測結果を事後に高精度で再現できたとしても，今後起こるであろう地震に対して同様のモデルを事前に構築することは困難であるから，工学的メリットは少ない．それに対して，比較的シンプルなモデルで観測結果を一定の精度で再現できるならば，同様のモデルを将来の地震に対して適用できる可能性が出てくる．本研究の目指すところはここにある．したがって，本研究で提案する複合型の震源モデルは，既往の研究（例えば田中他，2018；貴堂・永野，2020）でフリングス
テップの再現に用いられているモデルと比較して，よりシンプルなものとなっている．

図-2 複合型震源モデルの概念図

2. 複合型震源モデル
2.1 複合型震源モデルの全体像
本研究で提案する複合型震源モデルは，従来から設計において考慮されているアスペリティの破壊に伴う強震動に加え，フリングステップも考慮できるものである．図-2に複合型震源モデルの概念図を，表-1にモデル化の基

図-3 上盤側と下盤側の変位

表-1 複合型震源モデルの基本的考え方
モデル化の方法

主な計算対象領域B

一様な断層モデルに動的な情報を追加（主に浅部のすべりに起因する） 半無限地盤を対象とフリングステップ
した離散化波数法

領域A

深部にアスペリティを配置

深部のアスペリティに起因する
強震動

-5-

計算手法

修正経験的グリーン関数法（従来手法）

で計算するが，領域Bからの地震動は離散化波数法

ルス力を作用させた場合の波動場，すなわち半無限地盤（Bouchon, 1981; Bouchon, 2003）で計算し，その結果を加および水平成層地盤のグリーン関数を計算することがで算する．領域Bからの地震動を離散化波数法で計算するのきる．また，これをもとに，せん断食い違い型点震源によは，フリングステップの計算において，地表面における反る波動場，すなわち点震源による地震動を計算すること

射係数を厳密に扱うことが必要だからである．例えば図-ができる．さらに，点震源による地震動を多数重ね合わせ

3に示すような逆断層を考えた場合，例えば鉛直方向のズることで，面的な広がりをもった震源に対する地震動をレが1mとなる場合であっても，上盤側が0.5m上昇し，下計算することができる．このとき，3.で述べるように，断盤側が0.5m下降して，その和が1mとなるわけではない．

層面上に十分に密に点震源を配置することにより，永久自由表面の影響により上盤側の永久変位の方が大きいこ

変位成分を正確に含む地震動を求めることができる．

とが一般的である．その割合は，断層の傾斜角（図-3の𝛿𝛿）

離散化波数法では，震源から放射される地震波を様々やすべり角（断層面上で上盤側が下盤側に対して相対的な周波数と様々な波数𝑘𝑘を有する地震波の重ね合わせでにどちら向きに動くかを表す角度で，0度のときに純粋な表現する．ここで言う波数とは，厳密に言えば波数の水平

左横ずれ断層，90度のときに純粋な逆断層，180度のとき成分のことである．波数の小さな地震波はテイクオフ・アに純粋な右横ずれ断層，270度のときに純粋な正断層であングル（進行方向と水平面のなす角）の大きい波に対応しる）などに依存して決まってくる．一方，修正経験的グリており，波数の大きな地震波はテイクオフ・アングルの小ーン関数法などの一般的な強震動シミュレーション手法

さい波に対応している（図-4）．つまり，離散化波数法ででは，まず，震源断層を挟んで対称な地震動が生成され，は，震源から放射される地震波を様々な周波数と様々な

これが自由端反射条件により単純に2倍されると考える

テイクオフ・アングルを有する波に分けて考える．水平なことが一般的なので，上盤側と下盤側の永久変位の違い

地層境界に入射した地震波は透過もしくは反射するが，を表現することは困難である（サイト特性の違いは考慮

このとき，地層境界が水平であれば地震波の波数は変化できるが）．近地項や中間項の影響を考慮できる発展的なしない性質がある．そのため，半無限地盤および水平成層手法（野津，2006）を用いたとしても事情は同じである．

地盤では，地震波の透過と反射は周波数と波数に対して上盤側と下盤側の永久変位の違いを正しく計算するため独立に計算できる．そこで，周波数毎，波数毎に透過と反には，自由表面における境界条件を厳密に扱うことので射の計算を行い，その結果をすべての周波数と波数につきる手法を用いることが本質的に重要である．このこと

いて重ね合わせれば，任意の地点でのグリーン関数が求が，領域Bからの地震動を離散化波数法で計算する理由でまる．以上が離散化波数法の概要である．離散化波数法のある．

詳細については，最初の提案者であるBouchon自身の論文が参考になるほか，より詳しい定式化については野津本研究では，震源モデルをできるだけ簡便で利用しやす（2018b；2018c）も参考になる．

いものとするとの観点から，領域Bとしてはできるだけシンプルなものを考える．すなわち，個々の地震を表現する

さて，離散化波数法では，解が波数に関する無限級数のために複数枚の面を考えることはあっても，個々の面内形で与えられるため，実際の計算では無限級数をどこかではすべり量やすべり速度時間関数は一様と考える．まで打ち切らなければならない．その際の打ち切り波数

た，同様に，震源モデルをできるだけ簡便で利用しやすい𝑘𝑘𝑚𝑚𝑚𝑚𝑚𝑚 の決め方が，離散化波数法における誤差要因の一つものとするとの観点から，領域Bからの地震動は半無限地となり得る．フリングステップを正しく計算するために盤を対象とした離散化波数法で計算する．その際，領域Bは，震源が浅い場合の計算が重要であるが，震源が浅い場

からは0.5Hz以下の周波数成分のみが生成されると仮定し，合，打ち切り波数に対してより厳しい条件が求められる

0.5Hz以下の周波数成分のみを計算の対象とする．

ため，打ち切り波数に関する検討を慎重に行う必要があ

上述のような一様な断層モデルと半無限地盤の組み合る．また，面的な広がりをもつ断層は小断層の集合で表現

わせによりフリングステップを十分に再現できることはされるが，その際の小断層のサイズの決め方も，離散化波

熊本地震本震のデータから確認できる．この点について

数法における誤差要因の一つとなり得る．フリングステ

は4.で述べる．

ップを正しく計算するためには，断層面が近い場合の計算が重要であるが，断層面が近い場合，小断層サイズに対2.2 離散化波数法の概要

してより厳しい条件が求められるため，小断層サイズに離散化波数法（Bouchon, 1981; Bouchon, 2003）では，半関する検討を慎重に行う必要がある．

無限地盤および水平成層地盤の任意の位置に単位インパ-6-図-4 震源から放射される球面波（左）と様々なテイクオフ・アングルをもつ波の重ね合わせ（右）る知見を得ることを試みる．

以上のことから，3.においては，離散化波数法の誤差要因のうち，波数に関する和を計算する際の打ち切り波数

なお，離散化波数法では，震源から放射される地震波を𝑘𝑘𝑚𝑚𝑚𝑚𝑚𝑚 と小断層サイズの二点に着目し，離散化波数法によ様々な波数を有する地震波の重ね合わせで表現する際，る計算結果を静的な永久変位に関する解析解等と比較す

媒質が非減衰であるならば，特定の波数に対応する項のることで，フリングステップを正しく計算するために必

分母が0となってしまい，精度良く計算できなくなる問題

要な計算条件について検討する．

がある（野津，2018b）．そこで，離散化波数法が提案された当初（Bouchon, 1981）から，周波数に僅かに虚部を導入（Phinney, 1965）して分母が0となることを回避すること

3. 離散化波数法の計算条件が提案されていた．本研究の離散化波数法でも一貫してこの方法を採用している．その際の計算条件として，角周
波数𝜔𝜔に導入する虚部の大きさは𝜆𝜆 = 2𝜋𝜋⁄𝑇𝑇とした．ここにフリングステップを含む断層近傍地震動への離散化波

𝑇𝑇は離散フーリエ変換に用いるデータ長である．この方法

数法の適用性に関する検証（verification）の事例としては，離散化波数法による計算結果と永久変位に関する（静的を用いる場合，計算過程の最後にフーリエ逆変換で時刻な）解析解（Okada, 1985）との比較を行った事例（Honda

歴波形を得る際に波形の尾部にノイズが生じる（野津，and Yomogida, 2003）があるが，震源深さが100m未満とい2018s）．そこで，本研究では，あらかじめ必要なデータ長った厳しい条件での検討はこれまで十分に行われていなの2倍のデータ長で計算しておき，後半を棄却したものをい．また，永久変位に至る過程を含む変位波形について解計算結果とした（例えば計算結果として長さ8sの波形を得析解との比較によるverificationを行った事例は著者らのるためには，離散フーリエ変換に用いるタイムウインド

知る限り存在しない．これは，永久変位に至る過程を含む

ウは2倍の16sとし，最終的に得られる波形の後半部分を棄変位波形についての解析解が得られている場合が少ない却した）．

ためと考えられる．なお，ここでは，断層面上での面積分を解析的に実施したものを解析解と定義している．

3.1 永久変位に関する静的な解析解との比較

そこで，前者については，震源深さが100m未満といっ

まず，半無限地盤における永久変位に関する岡田によた厳しい条件を考え，離散化波数法による計算結果と永る解析解（Okada, 1985）との比較を行った．その際に着目

久変位に関する静的な解析解（Okada, 1985）との比較を行したのは，波数に関する和を計算する際の打ち切り波数った．また，後者については，上記の意味での数少ない解

𝑘𝑘𝑚𝑚𝑚𝑚𝑚𝑚 と小断層サイズである．

析解の一つである増田・引間による解析解（増田・引間，打ち切り波数𝑘𝑘𝑚𝑚𝑚𝑚𝑚𝑚 については次のような考察が可能で2002；野津，2017）を用い，断層近傍地震動の計算手法とある．離散化波数法においては，無限級数の各項には指数しての離散化波数法のverificationを行い，波数に関する和関数が現れるが，この指数関数は，フリングステップの計を計算する際の打ち切り波数𝑘𝑘𝑚𝑚𝑚𝑚𝑚𝑚 や小断層サイズに関す算において重要となる低周波側の極限では𝑒𝑒 −𝑘𝑘𝑧𝑧𝑠𝑠 のような

-7-

X (N)

IV

I

V
VI

II

100 m

III

90o

IV (50, -50)
V (50, 0)
Y (E) VI (50, 50)

pointsource
0m
10

Z

90o
point
sourcehalf-space

Y (E)
0m
10

Z

X (N)

case CX (N)case B
receiver coordinates(unit:m):I (0, -50)
II (0, 0)
III (0, 50)

100 m

receivers

case AX (N)case D

45o

45o
Y (E)

point source

0m
10

0m
10

point source
0m
10

Z

Y (E)

0m
10

Z

図-5 計算対象とした点震源．Nは北，Eは東，𝑧𝑧 = 0は地表面を表す．赤丸は点震源（point source）を表す．ローマ数字のⅠ～Ⅵは計算地点（receivers）を表す．

図-6 岡田による解析解との比較（点震源に対する結果），縦軸：変位，横軸：𝑘𝑘max = ηπ⁄𝑧𝑧s とした場合のηの値，破線：
岡田による解析解，マーカー：離散化波数法の結果-8-

図-6（つづき）

X (N)free surface

receivers
fault A:fault B:
fault C:
2.0km * 1.0km 10.0km * 5.0km 20.0km * 10.0km

45o

X (N)

Z

ul
t

Width

Fa

1.0 m

w
id
th

fault planerupture start point
LengthMaterial parameters:
α=5.6 km/sec, β=3.2 km/sec,
ρ=2.7 g/cm3,
QP= 500, QS=250

Y (E)

Z

図-7 面的な広がりを有する震源．Nは北，Eは東，𝑧𝑧 = 0は地表面（free surface）を表す．左の図で青い長方形は断層面（fault plane），赤丸は破壊開始点（rupture starting point）を表す．▲は計算地点（receivers）を表す．右の図は断面図である．右の図にあるように断層の大きさについては2km×1km（断層A），10km×5km（断層B），20km×10km（断層C）の3通りを検討した．

-9-

Okada-NS
DWN-125m
DWN-100m
DWN-50m
DWN-25m-2

20
10Displacement (10

m)

30

60

0
Okada-UD
DWN-125m
DWN-100m
DWN-50m
DWN-25m

40
20
0
0

1

2

3

4

Time (s)

5

6

7

40

Okada-NS
DWN-100m
DWN-50m
DWN-25m
DWN-20m
DWN-12.5mm)

10

Displacement (10

-2

30
20

8

0
-10
80

Okada-UD
DWN-100m
DWN-50m
DWN-25m
DWN-20m
DWN-12.5m

60
40
20
0
-20
0

1

2

3

4

Time (s)

5

6

7

40

Okada-NS
DWN-125m
DWN-100m
DWN- 50m
DWN- 25mm)

10

Displacement (10

-2

30
20

8

0
-10
80

Okada-UD
DWN-125m
DWN-100m
DWN- 50m
DWN- 25m

60
40
20
0
-20
0

2

4

6

8

Time (s)

10

12

14

16

図-8 岡田による解析解との比較（面的な広がりを有する震源に対する結果），上から順に断層A, B, Cに対応，最も断層に近い出力点での結果，縦軸：変位，横軸：時間，Okada：岡田による解析解，DWN：離散化波数法の結果
による結果である．この結果から，η = 5とすれば十分で形をしている．ここに𝑧𝑧𝑠𝑠 は震源と計算地点の𝑧𝑧座標の差のあること，すなわち，無限級数を打ち切る際の打ち切り波

絶対値である（計算地点が地表の場合は震源深さ）．そこで，打ち切り波数𝑘𝑘𝑚𝑚𝑚𝑚𝑚𝑚

としては，𝑒𝑒 −𝑘𝑘𝑧𝑧s が十分に小さくな

数は𝑘𝑘max = 5π⁄𝑧𝑧s とすれば十分であることがわかる．

次に，面的な震源を対象とする場合の小断層サイズにるような𝑘𝑘，すなわち𝑘𝑘𝑧𝑧s が十分に大きくなるような𝑘𝑘を選べば良いはずである．そこで，𝑘𝑘max = ηπ⁄𝑧𝑧s とし，このηついて検討を行った．図-7に示すように走向90°，傾斜45°，すべり角90°の長方形の震源を対象とし，図の▲のをどの程度の値とすれば良いか検討した．

図-5は検討に用いた点震源の位置と発震機構を示して位置（500m間隔）での地震動を離散化波数法で計算し，いる．図に示すように極めて浅い震源を対象とした．媒質

岡田による解析解との比較を行った．媒質のパラメターのパラメターは密度2.5×103kg/m3，P波速度5.0km/s，S波は密度2.7×103kg/m3，P波速度5.6km/s，S波速度3.2km/s，

速度3.0km/s， Qp=500，Qs=250である．すべりの時間関数

Qp=500，Qs=250である．断層の大きさについては2km×

はライズタイム1sのランプ関数とした．図のローマ数字の1km（断層A），10km×5km（断層B），20km×10km（断層ⅢとⅥの位置における各成分の永久変位に関して，離散C）の3通りを検討した．なお，無限級数を打ち切る際の打化波数法による結果と岡田の解析解を比較したものを図

ち切り波数は𝑘𝑘max = 5π⁄𝑧𝑧s とした．断層A, B, Cに対し，そ-6に示す．図の横軸は𝑘𝑘max = ηπ⁄𝑧𝑧s とした場合のηの値でれぞれ最も断層に近い出力点での結果を図-8に示す．

あり，破線が岡田による解析解，マーカーが離散化波数法

各々の断層面距離（断層面－出力点間の最短距離）は146m，- 10 -

25m，50mであった．図-8より永久変位を正確に含む地し破壊伝播の効果は入っていないこと，(2)すべり方向成震動を計算するためには，小断層サイズを断層面距離の分に関する解析解であり断層直交成分は対象外であるこ

半分以下とする必要のあることがわかる．

と，などの限界はあるものの，永久変位成分まで含んだ地震動の解析解として貴重な存在である．そこで，以下においては，解析解と同様の条件での地震動を離散化波数法

3.2 増田・引間による解析解との比較
永久変位の値のみならず，永久変位に至る過程を含むにより計算し，その結果を解析解と比較することにより，変位波形を対象とし，断層面上での面積分を解析的に実断層近傍地震動の計算手法としての離散化波数法の行している数少ない解析解として，限られた条件に対す

verificationを行う．離散化波数法は，本来，地表面でのるものではあるが，増田・引間による解析解（増田・引間，

地震波の反射や地層境界での地震波の透過・反射を考慮

2002；野津，2017）がある．これは，円形断層（図-9）でできる手法であるが，比較の対象とする解析解は全無限同時にすべりが生じたときの全無限弾性体における地震

弾性体に対するものであるため，離散化波数法の計算で動に関する解析解であり，円の中心軸に沿った観測点（図も震源からの入射波（野津，2018b；2018c）のみを考慮す-9）での変位のすべり方向成分を与えるものである．解析る．解析解の性質上，このような限られた条件での解は次式で与えられる．

verificationとならざるを得ないが，断層面の分割の細3 𝑎𝑎2 𝑟𝑟0 𝑅𝑅⁄𝛽𝛽⎡−� 𝜏𝜏 𝐷𝐷(𝑡𝑡 − 𝜏𝜏)𝑑𝑑𝑑𝑑 ⎤5⎢ 2 𝑅𝑅 𝑅𝑅⁄𝛼𝛼⎥⎢ 1 𝑎𝑎2 𝑟𝑟0 1⎥⎢−
𝐷𝐷(𝑡𝑡 − 𝑅𝑅 ⁄𝛼𝛼) ⎥32⎢ 2 𝑅𝑅2 𝛼𝛼⎥⎢+ 1 𝑎𝑎 𝑟𝑟0 1 𝐷𝐷(𝑡𝑡 − 𝑅𝑅 ⁄𝛽𝛽) ⎥⎢ 2 𝑅𝑅3 𝛽𝛽2
⎥
⎢ 1 1
⎥
⁄
)
𝐷𝐷(𝑡𝑡
−
𝑟𝑟
𝛽𝛽
⎢+
⎥
0
2
⎢ 2 𝛽𝛽⎥⎢ − 1 𝑟𝑟0 1 𝐷𝐷(𝑡𝑡 − 𝑅𝑅 ⁄𝛽𝛽)⎥
⎣ 2 𝑅𝑅 𝛽𝛽2⎦

かさや，波数に関する和を計算する際の打ち切り波数等に関する知見を得ることができると考えられる．
断層面と観測点の位置関係は図-9に示した通りである．
媒質の密度は2.7×103kg/m3，P波速度は5196m/s，S波速度(1)は3000m/sとした（P波速度はS波速度からポアソン固体を仮定して求めた）．円形断層の半径は1000mとした．すべりの時間関数としては指数関数型のものを用いた．

𝐷𝐷0𝐷𝐷(𝑡𝑡) = �1−𝑒𝑒 −1

ここに𝑟𝑟0 は観測点の𝑧𝑧座標，𝛼𝛼はP波速度，𝛽𝛽はS波速度，𝐷𝐷(𝑡𝑡)

はすべりの時間関数，𝑎𝑎は円形断層の半径，𝑅𝑅は円形断層

0�1 − 𝑒𝑒 −𝑡𝑡/𝑡𝑡𝑟𝑟 �𝐷𝐷0(𝑡𝑡 < 0)(0 ≤ 𝑡𝑡 ≤ 𝑡𝑡𝑟𝑟 )(𝑡𝑡𝑟𝑟 < 𝑡𝑡)(2)

ここに𝐷𝐷0 は最終すべり量であり2mとした．𝑡𝑡𝑟𝑟 はライズタの縁から観測点までの距離である．この解析解は，観測点イムであり1sとした．

を断層面に近づけていったときすべり時間関数の1/2に収面積分を数値的に実行するための点震源の座標は次式

束し（断層面を挟んで両側が変位するので1/2が付く），観
測点が十分に遠方の時は一般に期待されるように1⁄𝑅𝑅 にで与えた．

比例するという合理的な性質を有している．すべりの時間関数としては任意のものを与えることができる．
この解析解は，(1)円形断層が同時にすべった場合に対応

図-10

離散化波数法の計算に用いた点震源配置の例（点震源の間隔が 50m の場合）

図-9 対象とする円形断層（circular fault）と観測点（receiver）- 11 -

図-11

断層面距離が 10000m，1000m，100m の場合の変位波形のすべり方向成分の解析解と数値解の比較（𝑘𝑘max =5π⁄𝑧𝑧s ，点震源間隔は 50m），縦軸：変位，横軸：時間，黒線：増田・引間による解析解，赤線：離散化波数法の結果𝑥𝑥𝑖𝑖𝑖𝑖 = 𝑟𝑟𝑖𝑖 cos 𝜃𝜃𝑗𝑗
ここに𝑦𝑦𝑖𝑖𝑖𝑖 = 𝑟𝑟𝑖𝑖 sin 𝜃𝜃𝑗𝑗

𝑟𝑟𝑖𝑖 = (𝑎𝑎⁄𝑁𝑁𝑁𝑁 )𝑖𝑖𝜃𝜃𝑗𝑗 = (2𝜋𝜋⁄𝑁𝑁𝑁𝑁)𝑗𝑗

1
𝑁𝑁𝑁𝑁 = �6𝑖𝑖( 𝑖𝑖 = 0,1, … , 𝑁𝑁𝑁𝑁)( 𝑗𝑗 = 1,2, … , 𝑁𝑁𝑁𝑁)⋯ 𝑖𝑖 = 0⋯ 𝑖𝑖 ≥ 1(3)

である．𝑁𝑁𝑁𝑁は半径方向の分割数，𝑎𝑎⁄𝑁𝑁𝑁𝑁は点震源間隔であ(4)る．点震源の集合の例を図-10 に示す．ここでは点震源の間隔が 50m の場合を示している．このとき点震源の個数
は 1261 である．波数に関する和を計算する際の打ち切り

(5)
(6)

波数𝑘𝑘max は，3.1 の結果を踏まえ 𝑘𝑘max = 5π⁄𝑧𝑧s とした．ナ
イキスト周波数は 8Hz とした．

断層面距離が 10000m，1000m，100m の場合の変位波形(7)のすべり方向成分の解析解と数値解を比較して図-11 に- 12 -

示す．点震源間隔は 50m の場合を示している．10000m に

4.1 国土地理院による震源断層モデルとその適用おける変位波形は永久変位成分が少なく遠地項が卓越し断層運動による永久変位は，西原村小森（図-13）のよていると考えられるのに対し，100m における変位波形は

うに断層面距離のごく小さい地点では断層面全体のすべ

永久変位成分を顕著に含むものとなっている．10000m のりではなく観測点直近でのすべりによって決まるが，大計算結果にはナイキスト周波数が 8Hz であることに起因

津（KMM005）（図-13）のようにやや離れた場所では断層する誤差がやや見られるものの，いずれの地点において面全体の影響を受ける．ここでは，断層面からやや離れたも変位波形の計算精度は良好である．

場所でのフリングステップの再現も視野に入れ，断層面以上の結果を総合すると，離散化波数法によりフリンの浅部のみならず深部を含めたモデル化を行う．

グステップを正しく計算するためには，波数に関する和

断層運動による各地の永久変位をある程度説明でき，かを計算する際の打ち切り波数𝑘𝑘max は𝑘𝑘max = 5π⁄𝑧𝑧s あるいはつシンプルで扱いやすい震源モデルとして，国土地理院

それ以上，小断層サイズは断層面最短距離の1/2以下とすが地震直後に発表した震源断層モデル（国土地理院，2016）

ることが適切と考えられる．

がある．これは，図-13に示すように3枚の断層面からなるモデルである．これらは国土地理院の資料では西から順にB，A1，A2と名付けられているが，本稿では混乱を避け

4. 熊本地震本震への適用

るため断層1，断層2，断層3と呼ぶ．各断層の断層パラメ
ターを表-2に示す．ここでは，この震源断層モデルを複合型震源モデル構築の出発点とした．

次に，実地震である熊本地震本震を対象に，アスペリテ
ィの破壊に伴う強震動とフリングステップの両者を考慮

ここで，国土地理院による震源断層モデルを出発点とできる複合型の震源モデルを構築し，強震動シミュレー

することの利点はいくつか存在する．一つは，震源周辺のションを実施した．その際，震源モデルをできるだけ簡便

各地において，少なくとも永久変位をある程度適切に計で利用しやすいものとすることを念頭におき，測地学的算できることが保証されていることから，これに動的なデータ（静的な永久変位のデータ）を説明するように設定

情報（各断層面の相対的な破壊開始時刻，各断層面の破壊された国土地理院の震源モデルを出発点とし，これに最伝播様式，破壊伝播速度，ライズタイム）を付与すること

小限の動的な情報を付与することで領域Bを設定し，これで，各地のフリングステップを適切に計算できる可能性に領域Aを加えることで複合型の震源モデルを構築した．

があるという点である．このことは，断層モデル構築にお複合型震源モデルの構築手順を図-12に示す．

ける労力を大幅に低減してくれる．もう一つの利点は，国土地理院のモデルが，各断層面におけるすべり量を一様とした比較的シンプルなモデルであることから，工学的に扱いやすいという点である．仮に断層面上のすべり量の不均質性を考慮した複雑なモデルにより各地の永久変位とフリングステップを高精度で再現できたとしても，同様のモデルを今後起こるであろう地震に対して事前に構築することは困難であるから，工学的メリットは少ない．それに対して，比較的シンプルなモデルで各地の永久変位とフリングステップを一定の精度で再現できるならば，同様のモデルを将来の地震に対して適用できる可能性が出てくる．このような考えから，国土地理院による震源断層モデルをもとに，フリングステップを計算するための断層モデルを構築した．これを複合型震源モデルの領域Bとする．
ただし，後に見るように，すべりの一様な断層面を考え図-12るだけでは速度波形や加速度波形は過小評価となるので，複合型震源モデルの構築手順

これに領域A（アスペリティ）を加えることで，広い周波数帯域での地震動を計算できる複合型の震源モデルを構築した．

- 13 -

33.1★

Strong motion stations
Epicenter

Lat.(deg.)

33.0

KMM00532.9

南阿蘇村河陽
(KWY)
西原村小森
(KMR)

KMM00632.8

KMMH16

★

32.7

32.6

32.5
130.4

50km

0

130.6

130.8

131.0

131.2

Lon.(deg.)
図-13

熊本地震本震を対象とした複合型震源モデル．黒の長方形は領域 B であり，国土地理院による震源断層モデル（国土地理院，2016）の 3 枚の断層面に対応する．これらは，国土地理院の資料では西から順に B，A1，A2 と呼ばれているが，本稿では混乱を避けるため断層 1，断層 2，断層 3 と呼ぶ．赤の長方形は領域 A（アスペリティ）である．
▲は強震観測点（strong motion stations）を，★は震央（epicenter；破壊開始点を地表面に投影したもの）を表す．

表-2

国土地理院による震源断層モデル（国土地理院，2016）のパラメター（断層2のみ0.3kmだけ上方に移動して用いた）
MW
経度
緯度
上端深さ 長さ
幅
走向
傾斜
滑り角
滑り量
（km） （km） （°）（m）（°）（°） （km）（°）（°）
130.807 32.770 0.8
10.2
13.0
205
72
176
2.7
6.65
断層 1
130.996 32.878 0.6
20.0
12.5
235
60
209
4.1
6.96
断層 2
130.975 32.883 0.2
5.1
6.6
56
62
178
3.8
6.36
断層 3

以下の強震動シミュレーションでは，震源付近の媒質対象とした離散化波数法で計算する．ここで離散化波数のパラメターとして，密度2.7×103 kg/m3，S波速度3.5km/s，法を選択するのは，2.でも述べたように，逆断層や正断層P波速度6.1km/sを用いる．P波速度の6.1km/sはポアソン固の場合に下盤側よりも上盤側で永久変位が大きい性質を体を仮定してS波速度から求めた値である．これらの値か正しく考慮できるためである．また，計算精度を確保する

らせん断弾性係数𝜇𝜇を計算すると3.31×1010 N/m2となる．

ため，断層面距離（断層面－観測点間の最短距離）の半分

これと表-2に示す面積・すべり量から地震モーメントを

以下のサイズの小断層を用いる．また，領域Bからは0.5Hz

計算しモーメントマグニチュードに換算するとほぼ表-2以下の周波数成分のみが生成されると仮定し，0.5Hz以下に示すモーメントマグニチュードに対応する．

の成分を離散化波数法の計算対象とする．この点も2.で述べた通りである．熊本地震におけるフリングステップ

2.で述べたように，領域Bによる地震動は半無限地盤を- 14 -

200

200

KMMH16(NS)KMMH16(EW)
150

Displacement (cm)Displacement (cm)

150
100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)20

KMM006(NS)KMM006(EW)
150

Displacement (cm)Displacement (cm)

150
100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)20

KMM005(NS)150

Displacement (cm)150

Displacement (cm)

40

200

KMM005(EW)

100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)200200

KMR (EW)

KMR (NS)
150

Displacement (cm)150

Displacement (cm)30Time (s)

200

100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)5050

KWY (EW)

KWY (NS)0

Displacement (cm)0

Displacement (cm)

40

200

200

-50
-100
-150
-200

-50
-100
-150
-200

0

10

20

30

40

Time (s)図-1430Time (s)

0

10

20

30

40

Time (s)領域Bのみを考慮した変位波形の計算結果，縦軸：変位，横軸：時間，黒：観測結果，赤：計算結果（時間軸の原点
は気象庁の発震時刻）- 15 -

2と断層3の破壊開始時刻を，各々の破壊開始点と断層1の

50

NS-Displacement (cm)KWY

破壊開始点を結ぶ最短距離を破壊伝播速度で除した値とし，破壊伝播速度に関して試行錯誤を行ったところ，破壊0

伝播速度を2.4km/sとしたときに，西原村小森および南阿蘇村河陽におけるフリングステップの到来時刻が観測と-50

比較的整合したので，破壊伝播速度を2.4km/sに固定した．
その結果，断層2と断層3の破壊開始時刻はそれぞれ3.3s，-100

11.0sとなった．断層1の破壊開始時刻については，当初は0.0sとしていたが，この条件ではKMM006のNS成分のフリ-150

ングステップの到来が早すぎる結果となったので，2.5sと
した．

-200

すべり速度時間関数は次式のものに固定した．

-200

-150

-100

-50

0

50

𝐷𝐷

0
𝑉𝑉(𝑡𝑡) = (1−𝑒𝑒 −1
𝑒𝑒 −𝑡𝑡/𝑡𝑡𝑟𝑟 (0 ≤ 𝑡𝑡 ≤ 𝑡𝑡𝑟𝑟 ))𝑡𝑡

EW-Displacement (cm)図-15領域Bのみを考慮した変位軌跡の計算結果，縦軸：

𝑟𝑟

(8)

ここに𝐷𝐷0 は最終すべり量，𝑡𝑡𝑟𝑟 はライズタイムである．なお，

変位の南北成分，横軸：変位の東西成分，黒：観測結果，赤：計算結果

式(8)に示す関数は，入倉他（1997）に基づいて波形合成を行う場合で，かつ，分割数を十分に大きな値とする場合にの観測結果の分析から，0.5Hz以下の成分を計算すること

すべり速度時間関数が収束する関数である．この条件のでフリングステップは十分に再現されると判断した．

もとでライズタイムについて試行錯誤を行った結果，ラ

まず，国土地理院の断層モデルに対して，仮の動的なパ

イズタイムを2.0sとする場合に，観測波形を比較的再現す

ラメターを与え，図-13の各地点における変位波形を計算る結果となったので，ライズタイムとしては2.0sを採用ししたところ，最も断層に近い西原村小森では，EW成分のた．

永久変位が過小評価となった．そこで，断層2のみを上方

以上の条件の下で，領域Bのみを考慮して計算される各に移動させることを検討し，試行錯誤の結果，断層2を地の変位波形を観測結果と比較して図-14に示す．各地点

0.3kmだけ上方に移動させると，西原村小森でのEW成分におけるフリングステップの到来時刻は概ね適切に再現の永久変位が観測と整合したので，この状態で断層面のされている．したがって，破壊伝播速度の設定は概ね妥当

位置を固定した．なお，この結果として，後に示すように，であると考えられる．また，各地点において変位波形が立KMMH16のEW成分の永久変位が観測よりも大きくなっ

ち上がりに要している時間についても概ね再現されていたが，熊本地震の強震観測結果に占める西原村小森の記る．したがって，ライズタイムとして一様な値（2.0s）を録の重要性を考慮し，上方に移動させた断層面をこれ以用いたことは妥当であったと言える．これは，今後の地震降では用いる．表-2に示す静的なパラメターのうちで変におけるフリングステップを予測する上で好材料である

更を加えたのはこの点だけである．

と言える．なぜなら，各地の波形の立ち上がり時間を適切

次に動的なパラメターについて検討した．まず各断層に再現するために場所毎に異なるライズタイムを用いな面の破壊開始点について検討した．日奈久断層に相当すければならないとすれば，フリングステップの予測は困る断層1については，気象庁発表の震源から（正確には気難を極めると考えられるからである．

象庁発表の震源を通る鉛直線と断層1の交点から）破壊を

図-15には南阿蘇村河陽における変位軌跡の計算結果開始するものとした．断層2については，既往の波形インを観測結果と比較して示している．河陽ではまず西向きバージョン結果（例えばNozu and Nagasaka, 2017）を参考の変位が生じ，つづいて南向きの変位が生じているが，そにすると深部から破壊を開始したと考えられるため，北の状況が，断層2，断層3の2枚の断層面の動きにより再現

西側の端部から破壊を開始するとした．断層3については，できている．このように，測地学的データをもとに得られ最も震央に近い北西側の端部から破壊を開始するとした．

たシンプルな震源断層モデルに動的な情報を付与するこ各断層面では同心円状に破壊が進行するとした．

とで，ハンドリングの難しい不均質断層モデルを用いな

次に破壊伝播速度と各断層面の（震源での最初の破壊

くても，観測されたフリングステップをある程度再現でに対する）相対的な破壊開始時刻について検討した．断層

きることがわかった．

- 16 -

ところで，計算結果の速度波形を観測結果と比較して領域Bと領域Aを考慮して計算される各地の変位波形を

図-16に示すが，変位波形は概ね適切に再現できているに観測結果と比較して図-17に示す．変位波形の計算結果はも関わらず，速度波形は大幅な過小評価となっている．こ領域Bのみを考慮した場合に比べ永久変位の値に変化は

れは，ここでの計算ではアスペリティの寄与を考慮してないが少し短周期成分が乗っており観測波形の特徴に近いないためであると考えられる．そこで，4.2においては，づいている．領域Bと領域Aを考慮して計算される各地の領域Bに領域A（アスペリティ）を加えた複合型の震源モ

速度波形を観測結果と比較して図-18に示す．速度波形のデルを構築し，強震動シミュレーションを行う．

計算結果は領域Bのみを考慮した場合に比べ観測結果に大きく近づいている．
以上のように，実地震である熊本地震本震に対し，測地

4.2 アスペリティの考慮
著者らはすでに熊本地震本震を対象としたアスペリテ

学的データを説明するように設定された国土地理院の震ィモデルの作成を行っている（野津・長坂，2018）．ただ源モデルを出発点とし，これに最小限の動的な情報を付し，その際，アスペリティは布田川断層沿いの傾斜角84°

与することで領域Bを設定し，これに領域Aを加えることの断層面上に配置していた．それに対して，ここで領域B

で複合型の震源モデルを構築したところ，実際の観測事として設定している布田川断層沿いの断層2の傾斜角は実と整合するシミュレーション結果を得ることができた．

60°である（表-2）．そこで，ここでは傾斜角60°の断層
2上にアスペリティの再配置を行った．その結果設定され5. まとめたのが図-13に赤枠で示す領域Aである．設定した領域Aのパラメターを表-3に示す．領域Aからの地震動は修正経験的グリーン関数法（古和田他，1998；野津・菅野，2008；

本研究では，従来から設計において考慮されているア

野津他，2009）により計算し領域Bからの地震動に加算しスペリティの破壊に伴う強震動に加え，フリングステッ

た．領域Aからの地震動の計算手順の詳細については野プも考慮できる複合型の震源モデルを新たに考案した．

津・長坂（2018）を参照されたい．ここでは各地点におい考案した震源モデルは，アスペリティに対応する領域Aて経験的に得られているサイト増幅特性を用いている．

とその外側の領域Bからなる．領域Aは深部に留まるが，また中小地震記録の位相特性をグリーン関数の位相特性

領域Bは地表付近まで到達しフリングステップを生じる．

として用いている．西原村小森では2016年4月14日22:07に領域Aからの地震動は従来通り修正経験的グリーン関数

発生した地震の位相特性を用いており，南阿蘇村河陽で法で計算するが，領域Bからの地震動は離散化波数法で計は2016年4月16日1:45に発生した地震の位相特性を用いて算し，その結果を加算する．震源モデルをできるだけ簡便いる．

で利用しやすいものとするとの観点から，領域Bとしては，個々の地震を表現するために複数枚の面を考えることはあっても，個々の面内ではすべり量やすべり速度時間関表-3 領域Aのパラメター

破壊開始点
東経（deg）
破壊開始点
北緯（deg）
破壊開始点
深さ（km）
長さ（km）
×幅（km）
M0（Nm）
相対破壊開始時刻（s）
破壊伝播速度（km/s）
ライズタイ
ム（s）

アスペリテ
ィ1

アスペリテ
ィ2

アスペリテ
ィ3

130.817130.835130.872

32.83432.84432.855

法の計算条件については，離散化波数法による計算結果

9.8

9.8

8.1

久変位に至る過程を含む変位波形についての解析解（増3.0×3.0

3.0×3.0

5.0×5.0

0.20×1018

0.10×1018

2.0×1018

5.1

5.7

6.7

2.8

2.8

2.8

𝑘𝑘max = 5π⁄𝑧𝑧s あるいはそれ以上，小断層サイズは断層面最

0.27

0.27

1.00

用した．震源モデルをできるだけ簡便で利用しやすいも数は一様と考える．領域Bからの地震動の計算では半無限地盤を考え，0.5Hz以下の周波数成分のみを対象とする．
フリングステップを正しく計算するための離散化波数
を永久変位に関する静的な解析解（Okada, 1985）および永田・引間，2002）と比較することで確認した．離散化波数法によりフリングステップを正しく計算するためには，
波 数 に 関 す る 和 を 計 算 す る 際 の 打 ち 切 り 波 数 𝑘𝑘max は
短距離の1/2以下とすることが適切と考えられる．

複合型の震源モデルを実地震である熊本地震本震に適のとすることを念頭におき，測地学的データを説明するように設定された国土地理院の震源モデルを出発点とし，- 17 -

150

150

KMMH16(EW)

KMMH16(NS)100

Velocity (cm/s)

Velocity (cm/s)

100
50
0
-50

50
0
-50
-100

-100
-150

-150
0

10

20

30

0

40

10

Time (s)KMM006(EW)

50
0
-50

50
0
-50
-100

-100
-150

-150
0

10

20

30

0

40

10

Time (s)20

KMM005(NS)50

Velocity (cm/s)

50

Velocity (cm/s)

40

100

KMM005(EW)

0

-50

-100

0

-50

-100
0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)300

300

KMR (EW)

KMR (NS)

200

Velocity (cm/s)

200

Velocity (cm/s)

30

Time (s)

100

100
0
-100
-200

100
0
-100
-200

-300

-300
0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)150

150

KWY (EW)

KWY (NS)

100

Velocity (cm/s)

100

Velocity (cm/s)

40

KMM006(NS)100

Velocity (cm/s)

100

50
0
-50
-100

50
0
-50
-100

-150

-150
0

10

20

30

40

Time (s)

図-16

30

150

150

Velocity (cm/s)

20

Time (s)

0

10

20

30

40

Time (s)領域Bのみを考慮した速度波形の計算結果，縦軸：速度，横軸：時間，黒：観測結果，赤：計算結果（時間軸の原点
は気象庁の発震時刻）- 18 -

200

200

KMMH16(NS)KMMH16(EW)
150

Displacement (cm)Displacement (cm)

150
100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

20

200

KMM006(NS)150

Displacement (cm)Displacement (cm)

150
100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)20

KMM005(NS)150

Displacement (cm)150

Displacement (cm)

40

200

KMM005(EW)

100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)200200

KMR (EW)

KMR (NS)
150

Displacement (cm)150

Displacement (cm)30Time (s)

200

100
50
0
-50

100
50
0
-50

0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)5050

KWY (EW)

KWY (NS)0

Displacement (cm)0

Displacement (cm)

40

200

KMM006(EW)

-50
-100
-150
-200

-50
-100
-150
-200

0

10

20

30

0

40

Time (s)図-1730Time (s)Time (s)

10

20

30

40

Time (s)領域Bに加え領域Aを考慮した変位波形の計算結果，縦軸：変位，横軸：時間，黒：観測結果，緑：計算結果（時間軸の原点は気象庁の発震時刻）- 19 -

150

150

KMMH16(EW)

KMMH16(NS)100

Velocity (cm/s)

Velocity (cm/s)

100
50
0
-50
-100

50
0
-50
-100

-150

-150
0

10

20

30

0

40

10

Time (s)KMM006(EW)

50
0
-50

50
0
-50
-100

-100

-150

-150
0

10

20

30

0

40

10

20

100

KMM005(NS)50

Velocity (cm/s)

50

Velocity (cm/s)

40

100

KMM005(EW)

0

-50

-100

0

-50

-100
0

10

20

30

0

40

10

Time (s)

20

30

40

Time (s)300

300

KMR (EW)

KMR (NS)

200

Velocity (cm/s)

200

Velocity (cm/s)

30

Time (s)Time (s)

100
0
-100
-200

100
0
-100
-200

-300

-300
0

10

20

30

0

40

10

20

30

40

Time (s)Time (s)150

150

KWY (EW)

KWY (NS)

100

Velocity (cm/s)

100

Velocity (cm/s)

40

KMM006(NS)100

Velocity (cm/s)

100

50
0
-50
-100

50
0
-50
-100

-150

-150
0

10

20

30

0

40

Time (s)

図-18

30

150

150

Velocity (cm/s)

20

Time (s)

10

20

30

40

Time (s)領域Bに加え領域Aを考慮した速度波形の計算結果，縦軸：速度，横軸：時間，黒：観測結果，緑：計算結果（時間軸の原点は気象庁の発震時刻）- 20 -

ータを使用しました．記して謝意を表します．

これに最小限の動的な情報を付与することで領域Bを設（2022年2月1日受付）

定し，これに領域Aを加えることで複合型の震源モデルを構築したところ，実際の観測事実と整合するシミュレー
ション結果を得ることができた．特にフリングステップ

参考文献については，ハンドリングの難しい不均質断層モデルを

安積恭子・岡上政史・杉山裕樹・井上直人・北田奈緒子

用いなくても，領域Bの寄与としてある程度説明できるこ（2021）：大阪湾岸道路西伸部の橋梁設計において考とがわかった．

慮するとう曲変位量の検討，土木学会第76回年次学術講演会，I-366．
入倉孝次郎（1996）：阪神・淡路大震災－防災研究への取6. おわりにり組み－阪神淡路大震災を引き起こした強震動，京都大学防災研究所年報，No.39A，pp.17-33.

本研究では，従来から設計において考慮されているア

入倉孝次郎・香川敬生・関口春子（1997）：経験的グリースペリティの破壊に伴う強震動に加え，フリングステッ

ン関数を用いた強震動予測方法の改良，日本地震学

プも考慮できる複合型の震源モデルを新たに考案した．

会講演予稿集，No.2，B25.

特に，複数枚の一様断層モデルによりフリングステップ

貴堂峻至・永野正行（2020）：2016年熊本地震における大が説明可能であることが示された点は重要な結果である

振幅速度パルスの解釈と断層近傍地震動に対する浅と考えられる．仮に断層面上のすべりの不均質性を考慮

部・深部断層破壊の影響，日本建築学会構造系論文集，

した複雑な震源モデルを導入しない限りフリングステッ

No.773，pp.879-889．

プを事後に説明できないとするならば，そのような情報

国土地理院：平成28年熊本地震の震源断層モデル（暫定），が事前に得られることは期待できないから，フリングス

2016，http://www.gsi.go.jp/common/000140781.pdf

テップを予測し設計で考慮することは難しくなる．実際

古和田明・田居優・岩崎好規・入倉孝次郎（1998）：経験には，比較的シンプルな震源モデルにより，各地の永久変的サイト増幅・位相特性を用いた水平動および上下

位とフリングステップを一定の精度で再現できることが動の強震動評価，日本建築学会構造系論文集，No.514，

本研究でわかったので，今後は同様のモデルを予測や設

pp.97-104.
田中信也・金田惇平・引間和人・久田嘉章（2018）：地表計に活用できる可能性がある．
ただし，どのような構造物の設計においてフリングス

地震断層近傍における永久変位を含む長周期成分のテップを考慮すべきかという点については，少なくとも

地震動評価のための震源モデルの設定方法，日本建港湾分野ではこれまでほとんど議論されておらず，今後

築学会構造系論文集，No.752，pp.1525-1535．

の議論に待つ必要がある．また，フリングステップに限ら

豊岡亮洋・小野寺周・坂井公俊（2021）：慣性力と断層変ず，地震動の振幅の予測結果には当然ではあるが大きな位を同時に受けるラーメン高架橋の地震応答に関する基礎的検討，土木学会第76回年次学術講演会，I-364．

不確実性がある．この点については，従来の強震動シミュ
レーションでは，まずは断層面上でのすべり量を想定し，

野津厚（2006）：統計的グリーン関数法に近地項と中間項

強震動シミュレーションを行った結果として地震動の振を導入するための簡便な方法，第12回日本地震工学

幅が決まるという場合が多かった．しかしながら，フリン

シンポジウム論文集，pp.190-193．

グステップは地表における断層変位量と直接関係してお野津厚・菅野高弘（2008）：経験的サイト増幅・位相特性

り，断層変位量については，事業によっては設計値があらを考慮した強震動評価手法－因果性と多重非線形効

かじめ与えられている場合がある（例えば安積他，2021）．

果に着目した改良－，港湾空港技術研究所資料，

このようなことを考慮すると，従来の手順とは逆に，計算No.1173．

結果のフリングステップが所要の断層変位量と整合する

野津厚・長坂陽介（2018）：2016年熊本地震の前震と本震

ように断層面上でのすべり量を決めるという方法も考えの震源過程の推定と特性化震源モデル，港湾空港技られる．このような点についても今後十分に議論する必術研究所報告，Vol.57，No.1，pp.3-35.

要がある．

野津厚・長尾毅・山田雅行（2009）：経験的サイト増幅・位相特性を考慮した強震動評価手法の改良－因果性を満足する地震波の生成－，土木学会論文集A，Vol.65，

謝辞

pp.808-813.

防災科学技術研究所，気象庁，熊本県による強震波形デ- 21 -

野津厚（2017）：運動学的断層モデル（その1）（増田・引T., Miyashita, Y., Mori, H., Imanishi, K., Takeda, N., Ochi,

間 ， 2002 ）， https://www.pari.go.jp/bsh/jbn-kzo/jbn-

T., Otsubo, M., Asahina, D. and Miyakawa, A. (2016):

bsi/taisin/tutorial_jpn/tutorial_003.pdf

Characteristics of the surface ruptures associated with the2016 Kumamoto earthquake sequence, central Kyushu,

野津厚（2018a）：Phinney法，https://www.pari.go.jp/bsh/jbn-

Japan, Earth Planets Space, Vol.68, Article No.191.

kzo/jbn-bsi/taisin/tutorial_jpn/tutorial_014.pdf
野 津 厚 （ 2018b ）： 水 平 成 層 地 盤 の グ リ ー ン 関 数 ，
https://www.pari.go.jp/bsh/jbn-kzo/jbnbsi/taisin/tutorial_jpn/tutorial_020.pdf
野津厚（2018c）：水平成層地盤におけるせん断食い違い型
点震源による地震動，https://www.pari.go.jp/bsh/jbnkzo/jbn-bsi/taisin/tutorial_jpn/tutorial_023.pdf
林康裕・杉野未奈（2018）：パルス性地震動に対する建物の応答と地震対策，第15回日本地震工学シンポジウ
ム，pp.1025-1034.
増田徹・引間和人（2002）：断層近傍における地震動，日本地震学会講演予稿集，A16．
Bouchon, M. (1981): A simple method to calculate Green's
functions for elastic layered media, Bulletin of the
Seismological Society of America, Vol.71, pp.959-971.
Bouchon, M. (2003): A review of the discrete wavenumber
method, Pure and applied geophysics, Vol.160, pp.446-465.
Hisada, Y. and Bielak, J. (2003): A theoretical method forcomputing near-fault strong motions in layered half-spaceconsidering static offset due to surface faulting, with aphysical interpretation of fling step and rupture directivity,
Bull. Seism. Soc. Am., Vol.93, No.3, pp.1154-1168.
Honda, R. and Yomogida, K. (2003): Static and dynamic
displacement near a fault with the discrete wavenumber
method, Phys. Earth Planet Inter., Vol.137, pp.107-127.
Nozu, A. and Nagasaka, Y. (2017): Rupture process of the main
shock of the 2016 Kumamoto earthquake with specialreference to damaging ground motions: waveform
inversion with empirical Green's functions, Earth Planets
Space, 69:22, DOI 10.1186/s40623-017-0609-3.
Okada, Y. (1985): Surface deformation due to shear and tensile
faults in a half-space, Bull. Seism. Soc. Am., Vol.75,
pp.1135-1154.
Phinney, R.A. (1965): Theoretical Calculation of the Spectrumof First Arrivals in Layered Elastic Mediums, Journal ofGeophysical Research, Vol.70, pp.5107-5123.
Sekiguchi, H., Irikura, K. and Iwata, T. (2000): Fault geometryat the rupture termination of the 1995 Hyogo-ken Nanbu
earthquake, Bull. Seism. Soc. Am., Vol.90, No.1, pp.117133.
Shirahama, Y., Yoshimi, M., Awata, Y., Maruyama, T., Azuma,- 22 -

港湾空港技術研究所資料

No.1401２０２２．３

編集兼発行人

国立研究開発法人海上・港湾・航空技術研究所

発

港 湾 空 港 技 術 研 究 所

行

所

横 須 賀 市 長 瀬 ３ 丁 目 １ 番 １ 号
TEL. 046(844)5040 URL. http://www.pari.go.jp/

C （2021）by MPAT
Copyright ○All rights reserved. No part of this book must be reproduced by any means without the written
permission of the President of MPAT
この資料は、海上・港湾・航空技術研究所理事長の承認を得て刊行したものである。したがって、
本報告書の全部または一部の転載、複写は海上・港湾・航空技術研究所理事長の文書による承認を得ずしてこれを行ってはならない。