「強震動予測－その基礎と応用」第８回講習会

波数積分法による強震動計算

工学院大学・建築学科
久田嘉章

概要○強震動予測の必要性
標準地震波（告示波など）と様々な強震観測記録の違い（指向性パルス、フリングステップ）○波数積分法による強震動計算法
震源・伝播・サイト特性○計算コードと事例の紹介

特徴的な観測強震動(変位波形)と地震被害
2003年十勝沖地
震の苫小牧市におけるナフサタンクの全面火
災

堆積層表面波（長周期地震動）

震源近傍の強震動（指向性パルス）

1995年兵庫県南部地震における神戸市の木造家屋の倒壊

地表断層近傍の強震動（フリングステップ）

1999年台湾・集集地震による地
表断層上のRC
建物の傾斜

震源特性：東西圧縮場と断層タイプ
北東－南西 逆 断 層の走向
東傾斜
右横ずれ
北断層
上盤
下盤
東
（平野側） （山側）

東西
圧縮

右横ずれ断層南北走向
上盤
下盤（山側） （平野側）

東西
圧縮

左横ずれ断層北西－南東の走向

日本列島で支配的な東西圧縮場における断層走向と断層タイプ

逆断層

1995年阪神淡路大震災（M7.3）
Ｍ７、活断層地震、震度７

地震調査研究推進本
部：｢日本の地震活動｣

阪神大震災写真20選（左上） http://web.kyoto-inet.or.jp/org/gakugei/kobe/photo20.htm#01phot_h• 被害概要（平成１２年消防庁災害対策本部）
・死者：6,432名（約９割は建物倒壊）
・負傷者：４３,７９２名 ・建物被害：５１２,８８２棟

1940年Imperial Valley 地震（横ずれ断層）

El Centro波→ランダム波

1979年Imperial Valley 地震（横ずれ断層）

Meloland波→指向性パルス波

１９９５年兵庫県南部地震（横ずれ断層）
神戸大学（KBU）波→指向性パルス

１９９５年兵庫県南部地震
アスペリティーと指向性パルス

淡路側

神戸側

釜江・入倉（1997）

１９９５年兵庫県南部地震
アスペリティーと指向性パルス
JMA-simu(FN)
JMA-obs

V e l(k in e )

155
0-155
0510Time(sec)

15V e l(k in e )

155

2025

MOT-simu(FN)
MOT-obs

0-155
05Time(sec)

15No.5の寄与

V e l(k in e )155

松島・ 川瀬（2000）

102025
TKTsimu(FN)

( )

-220
0

510Time(sec)152025

1995年兵庫県南部地震における地震動（速度波形）の卓越方向

１９９５年兵庫県南部地震
木造家屋の倒壊方向

１９９５年兵庫県南部地震・神戸市で観測された強震記録の地震応答スペクトル
加速度応答スペクトル（5%）

3000

改正基準法（工学的基盤）
改正基準法（第一種地盤）改正基準法（第二種地盤）改正基準法（第三種地盤）
JMA神戸(FN)
サイト特性有り
JMA鷹取(FN)JMA葺合(FN)（1.5倍程度）

加速度応答（gal）

2500
2000
1500
1000
500
00

1

23

45周期（秒）逆断層と指向性パルス
パルスが パルスが パルスが出易い
出難い
出難い

パルスは出難い

震央

告示凌駕

断層面
直交成分

傾斜角小
傾斜角大

震源

震源

高角逆断層（例：1994年ノースリッジ地震、2004年新潟県中越地震）

低角逆断層（例：東京湾北西部地震）

1994年米国・Northridge地震（逆断層：震源モデル）
A’
断層面延長線

断層面
延長線

断層面＋15°線

5 km
断層面＋15°線

15°
15 km

A

震源

傾斜角

FN

40°

FP
Wald ほか(1996)に震源モデル

1994年Northridge地震（強震動）1000

Forward (FN)

1000

U56

800

SYL

400

SCS

200

velocity (cm/s)

NHL

600

NHL

600

SYL

400

SCS

200

RRS

001020RRS

03001020-200time (s)-200time (s)1000

Backward (FN)

1000

Backward (FP)

SSU

800

U55

600

U53400

U03

200

SVA

00-2001020time (s)

U55

600

U53400

U03200

SVA

0030
-200

30

SSU

800
velocity (cm/s)

velocity (cm/s)

800

velocity (cm/s)

Forward (FP)U56

1020time (s)301994年Northridge地震（速度応答）
Forward (FN)

U56
NHL
SYL
SCS
RRS

300
250
200

velocity response (cm/s

Velocity Response (cm/s

350

150
100
50350

Forward (FP)300

震源・地盤？

250
200
150
100
50

000246810024350

6810period (s)period (s)Backward (FP)350

Backward (FN)
velocity response (cm/s

Velocity Response (cm/s

U56
NHL
SYL
SCS
RRS

SSU
U55
U53
U03
SVA

300
250
200
150
100
50300

SSU
U55
U53
U03
SVA

250
200
150
100
50
0002468010246810

period (s)period (s)東京湾北西部地震(M7.3)の強震動計算
1500

1994ノースリッジ地震タイプ
Waldモデル(N26E)

長さ63.6km

北

東

観測点
（新宿）
震源301000

震源1

震源5

幅 31.2km

20km

加速度(gal)

震源2500

内閣府モデル(N26E)
T

0102030

40-500
断層面の投影図及び震源観測点配置図-1000

断層直交成分（N26E）
断層傾斜角=25度60Time(s)

震源4
図１50

JMA神戸

h=0.05＋σ
100波平均断層平行成分（N116E）－σ→ 断層モデル2通り×震源5通り×応力降下量分布2通り×ランダム数5通りの計100波形の強震動計算

１９９９年台湾集集地震と地表断層出現• 巨大な地表断層の出現→フリングステップ

日本大学理工学部のホームページより

１９９９年台湾集集地震と観測強震動

地表断層出現とフリングステップ
1999 台湾・集集地震（逆断層・上盤側）

1999 台湾・集集地震（地表断層近傍の地震被害・全戸調査）

1999 台湾・集集地震（地表断層近傍の地震被害・全戸調査）

1999 台湾・集集地震（地表断層近傍の地震被害・全戸調査）

海洋型地震と内陸型地震
糸魚川-静岡構造線
上町断層帯

すべり
10 m→3m

宮城県沖地震南関東地震
南海・東南海・東海地震

神縄・国府津－松田断層帯

富士川河口断層帯
すべり 7m

日本周辺の被害地震（｢日本の地震活動｣より）

全国を概観した地震動予測値図：地震危険度マップ（1 kmサイズ）
・プレート境界地震、活断層などの震源の発生確率、地盤増幅の影響を考慮
・震度（地震動の強さ）・考慮期間・確率の関係情報を表示
・地震波形もダウン
ロード可能

全ての地震を考慮し、30年以内に震度６弱以上の地震動が発生する確率

強震動予測と入力地震動：地震動予測値図の公開（例：JSHS）
350
300
250

250

200
速度応答（cm/s）

200
速度（cm/s）

糸井川－静岡構造線（北・中部断層）

NS

150
100
50

EW

150
100

改正基準法（安全限界）
糸井川-静岡(NS)500

糸井川-静岡(EW)

0-50

10

20

30

40

50

60

時間（秒）

00

2

4
周期（秒）

-100

主要９８断層の強震動予測波形

糸井川ー静岡構造線（北・中部）の速度波形例

速度応答スペクトル

6震源近傍における強震動（まとめ）１）ランダム位相波（エルセントロ波・告示波）
告示波レベル（短周期が卓越、継続時間はほぼ震源の破壊時間）→免震が有効
２）指向性パルス（神戸波・ノースリッジ波）
M7以下では通常、断層破壊の進行方向、断層面の直交成分に現れる破壊力が非常に大きいキラーパルス。震源アスペリティーより発生→M7以
下では免震は有効、M７以上は？
３）フリングステップ（台湾波・ランダース波）
地表断層のごく近傍における長周期で卓越する断層すべり→耐震が有効。変形対策が必要

波数積分法による強震動計算• 強震動地震学における理論的手法（文献１）• 表現定理とグリーン関数（文献１、２）• 波数積分法：平行成層地盤におけるグリーン関数（文献２）→高振動数でも計算可能→統計的震源モデルに適用（文献３）• 震源特性：指向性効果、地表断層効果（文献１）• 伝播・サイト特性：実体波と表面波の伝播・増幅（文献１）• 計算事例の紹介（文献１）
計算コード・データのダウンロード
http://kouzou.cc.kogakuin.ac.jp/ → Open Dataへクリック
参考文献
１） 久田嘉章、地盤震動－現象と理論－、第４章 成層地盤の解析手法（分担）、日本建築学会、84-102, 2005
２） 久田嘉章、建築の振動：応用編、６章 地震と地震動（分担）、80-140、朝倉書店、2008
３） Hisada, Y., Broadband strong motion simulation in layered half-space using stochastic
Green‘s function technique, Journal of Seismology, Vol.12, No.2, 265-279、2008

フーリエ変換
強震動地震学における定式は振動数領域で行なう
フーリエ変換∞F (ω ) = ∫ f (t )e+ iωt

−∞

dt

1 ∞
− iω tフーリエ逆変換 f (t ) =F (ω )e dω∫
2π −∞注意：工学におけるフーリエ変換の定義とは指数部の符号が逆！

強震動計算 → ハイブリッド手法
低振動数：理論的手法→波数積分法、差分法など
高振動数：経験的・統計的手法→統計的グリーン関数法など
マッチング振動数（0.5～2.0 Hz）

フリングステップ、指向性パルス、表面波・・・

ランダム波形・・・

波数積分法で扱える震源・地盤モデル
フリングステップ

表面波（ラブ波、レイリー波）

観測点
表層地盤

地表断層

工学的基盤
地中断層

地震基盤

実体波
（P波、S波）平行成層地盤モホ反射波
震源

モホロビッチ面

強震動計算の基礎式（表現定理）
表現定理（面震源による強震動の計算式）

U k ( Y ; ω ) = ∫ T ik ( X , Y ; ω ) [D i ( X ; ω ) ]d ΣΣ

観測点Yの変位解

グリーン関数 ソース点Xのダブルカッ
プル震源・すべり関数

ソース点 X

観測点 Y

グリーン関数
破壊フロント （地震動の伝播特性）
破壊開始点（震源）震源特性：点震源（ダブルカップル震源）

点震源（ラディエーションパターン）面震源モデルと震源パラメータ
地表面

X (北)
走向角

Y( 観測点)(ϑ )Y (東 )

断 層走向グリー ン関数（ 長さ方向）Z (下)ア スペ リ テ ィ

す べ り （D、上 盤で定義）

断層幅（W）
断層面の法線ベクトル（n）

すべり角
破壊伝播r

震源（破壊開始点）

断層長さ（L）

λX(ソース点)
破壊フロン トア ス ペ リテ ィ

下盤
傾斜角 (δ )断層走向（長さ方向）

擬似動力学的すべり速度関数

震源特性：震源近傍の強震動
ランダム波と指向性パルスの成因• 横ずれ断層
北

パルスが出難い

パルスが出易い

震源特性の計算例
地中断層における指向性効果

地盤は半無限一様地盤
断層面を５×２＝１０の小断層に分割
各小断層で３×３＝９点のガウス面積分を実施（滑らかな破壊伝播を再現）

震源特性の計算例
地中断層における指向性効果
1515

断層面平行成分（EW）

断層面平行成分（EW）

1010

断層面直交成分（NS）
0-50

10

20

30
時間（秒）

-10-155
速度（cm/s）

速度（cm/s）

5断層面直交成分（NS）

0-50

102030時間（秒）-10
観測点１-20-15

観測点２-20

指向性パルス
速度波形（観測点２）

速度波形（観測点１）

震源特性：震源近傍の強震動
地表断層とフリングステップ• 横ずれ断層
N

観測点１

観測点２
A

破壊フロントが遠ざかる

観測点３B

破壊フロントが近づく

断層面ごく近傍における表示定理の断層面積分の特異点の扱い
U k ( Y ; ω ) = ∫ T ikD ( X , Y ; ω ) [D i ( X ; ω ) ]d ΣΣ141210In te g ran d

86420

1

224549837

Length1541131721252933

95• 観測点が断層面近傍にある
時、グリーン関数が1/r2に近いオーダーで発散• 発散点近傍に多数のグリーン関数を分布させる必要• 動的な成層地盤では多大な計算時間• 特異性は静的項に起因

Width

1

グリーン関数の断層面分布→ 静的項で発散

断層面近傍における強震動計算（表示定理）
表示定理（オリジナル）

U k (Y ; ω ) = ∫ T ik ( X , Y ; ω ) [D i ( X ; ω ) ]d ΣΣ

修正表示定理（Hisada and Bielak, 2002）{}

U k (Y ; ω ) = ∫ Tik ( X , Y ; ω ) − TikS ( X , Y ) [Di ( X ; ω ) ]dΣΣ+ ∫ T ikS ( X , Y ) [D i ( X ; ω ) ]d ΣΣTik ：動的グリーン関数 → 幾何減衰: 1/r～1/√r → パルス波
TikS ：静的グリーン関数 → 幾何減衰: 1/r２ → Fling Step

震源特性の計算例
地表断層によるフリングステップ

震源特性の計算例（フリングステップ）
指向性パルス

フリングステップ

フリングステップ

強震動計算の基礎式（表現定理）
表現定理（面震源による強震動の計算式）

U k ( Y ; ω ) = ∫ T ik ( X , Y ; ω ) [D i ( X ; ω ) ]d ΣΣ

観測点Yの変位解

グリーン関数 ソース点Xのダブルカッ
プル震源・すべり関数

ソース点 X

観測点 Y

グリーン関数
破壊フロント （地震動の伝播特性）
破壊開始点（震源）

平行成層地盤のグリーン関数（伝播・サイト特性）○伝達マトリックス法
・古典：Haskellのマトリックス法→高振動数で発散
・現在：R/Tマトリックス法（Kennett, Lucoand Apselなど）→高振動数で安定

平行成層地盤のグリーン関数（波数積分法：フーリエ・ベッセル変換より）

ベッセル関数

震源の加振向き
観測点の振動向き

深さ方向の振幅 （SH波）（P－SV波）

波数=ω/C

波数軸上で関数値が発散する極が表面波、分岐点は実体波に対応

実体波
速度Vp

速度Vs表面波速度
位相速度（C）→波形の山谷Rayleigh wave

群速度（U）→波群（エネルギー）の速度
常にC≧ULove Wave

実体波と表面波の計算（２層地盤の例）

実体波と重複反射波（震源１と観測点１）

実体波と表面波の計算（２層地盤の例）

SH波とラブ波の伝播（震源２と観測点２～６）X成分

表面波：距離とともに波形がバラバラになる→分散

P・SV波とレイリー波の伝播（震源２と観測点２～６）

表面波：距離とともに波形がバラバラになる→分散分散曲線（モードと周期と位相速度・群速度）

単純な２層地盤（基盤と表層）の例

表面波の卓越振動数

ラブ波のミディアムレスポンス

レイリー波のミディアムレスポンス

波数積分法による強震動計算
公開計算コード（Fortran）○主な特徴• 伝達マトリックス：R/Tマトリックス（Luco & Apsel）→高振動数まで安定して計算• 表面波の位相速度（極）を事前に計算（群速度・分散曲線・変位応力ベクトル・H/Vスペクトル・Medium Responseも計算・・）• 振動数依存のQ値を導入可能• すべり速度関数：二等辺三角形の重ね合せでどんな形状にも対応• 断層面積分：小断層でガウス積分（1×1～6×6点）• 地表断層によるフリングステップ（静的項）も計算• 破壊フロントにランダムな乱れを導入可能

計算事例紹介：1992年Landers地震
Wald and Heatonのモデル(1994)
文献：Hisada, Y, and J. Bielak, Effects of Sedimentary Layers
on Directivity Pulse and Fling Step, 13WCEE, 2004

観測記録と計算波形の比較（LUC）
断層面直交

断層面平行

速度

速度

断層面直交
変位
断層面平行
変位

波数積分法による強震動計算の流れ○計算手順①表面波の位相速度（分散曲線）を計算（波数軸上の極に相当）：phs3sQ-v3.f②波数積分法による強震動計算（振動数領域）：grflt12sx1-v3.f③フーリエ逆変換（時間領域）：grfftsp.f

入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
時間・周波数データ： phs3sQ-v3.f、grflt12sx1-v3.f*** Data for Delta Time, Duration, and Minimum Period ***Delta TimeNumber of Time (must be Power of 2)
0.5
128
Minimum PImaginary Omega for Phinney's method
10.01→時間刻みを0.5秒で128点（=27、継続時間=64秒）で、最小周期１秒（＞時間刻みの２倍）で計算。円振動数ωの虚数部を0.01とする（Phinney法のため）。

地盤モデル（平行成層地盤）
Num. Density (g/cm3) Vp (km/s) Qp Vs (km/s) Qs Thickness (km)12.3

3.8100

1.98

30

1.52

2.65.5600

3.15300

2.53

2.7

6.26003.52300

22.04

2.87

6.86003.83300

6.05

3.5

8.06004.64300

0.0

入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
地盤データ： phs3sQ-v3.f、grflt12sx1-v3.f*** MEDIUM DATA ***NL (N (Note: Frequency-Dependent Q; Qp(f) = Qp0 x f ** Qpf & Qs(f) = Qs0 x f ** Qsf
5
Layer density(t/m3)Vp(m/s) Qp0 (Qp(f)=Qpf Vs(m/s) Qs0 (Qs(f) Qsf Thichness(m)12.33800
10001980
3001500
22.65500
600
03150
300
02500
32.76200
600
03520
300022000
42.876800
60003830
300
06000
53.58000
600
04640
300
00→５層地盤。Q値の振動数依存なし。
密度、Vp、Qｐ、Vs、Qｓ、層厚震源モデル：３つの断層面の重合わせCamp Rock/Emerson
Fault (CEF)

Homestead Valley
Fault (HVF)
Johnson Valley
Fault (JVF)

観測点(Lucerne Valley)

震源断層基準点・小断層番号ほか

入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
断層モデルデータ１： grflt12sx1-v3.f*** Seimic Fault Parameters (after Wald and Heaton, BSSA, 1994 Combined
Length (m) Width (m) Num. of SNum. of SNumber of Gaussian Points per Su
36000
15000
1264Start Time oStrike (deg) Dip (deg) Vr (m/s) dtr (s; av random number for dtr (in0320
90
2700
01
Location of Y(m)Z(m)(Note: X->North, Y->East, Z->Down)
21818
5119
7000
Location of Y(m)Z(m)(Note: X->North, Y->East, Z->Down)
50949
4015
15000

→断層長さ・幅、小断層分割数、ガウス点数
走向角、傾斜角、破壊伝播速度、フロント乱れ無震源位置、断層基準点位置(X,Y,Z)

すべり速度関数（小断層ごと）
二等辺三角形の重合せ（６個）すべり量（面積）
すべり速度

時間（秒）
破壊開始時間１秒

入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
断層モデルデータ２： grflt12sx1-v3.f
Number of TInterval Tim Slip Velocity Func.(
611Time Windo 1st Half Rise2nd Half Rise Time
1
0.5
0.5
2
0.5
0.5
3
0.5
0.5 すべり速度
4
0.5
0.5
5
0.5
0.5
6
0.5
0.5

すべり量（面積）→すべり速度関数
時間（秒）破壊開始時間１秒入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
断層モデルデータ３： grflt12sx1-v3.f
Slip and Rake Data Pattern (=1: Regular Patten (ex. grflt12s.f), =2: Time Window1Sub-FaultSlip of 1st Slip of 2nd Slip of 3rd Slip of 4th Slip of 5th Slip of 6th10.02
0.13
0.24
0.11
00.15
200.340.08000300.29000040.33
0.42
0.13
00
0.02
Sub-FaultRake of 1stRake of 2stRake of 3stRake of 4stRake of 5stRake of 6st
11801801801801801802180180180180180
180
3180180
180180180
180
4180180
180180180
180→小断層・すべり時間ごとのすべり量・すべり角

入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
波数積分データ： grflt12sx1-v3.f*** Data for Static Wavenumber Integration using Greenfield's Quadrature **The first cInitial Number of Intgegration Points for Adaptive Newton-Cotes Qua
216
The seconInitial Number of Intgegration Points for Adaptive Newton-Cotes Qua
1032*** Data for Dynamic Wavenumber Integration using Simpson's and Filon's quNumber of Integration Points from 0 to om/Ryleigh(min)200
Number of Integration Points from om/Ryl(min) to om/c(final)
50
Factor for c(final): c(final)=Ryl(min)/Factor (Ex., 10)10→静的項：積分路変換法＋Newton-Cotes則
動的項：Simpson則（小波数）、Filon積分法（大波数）入力データ例：grflt12s1.in-CEF（Camp Rock/Emerson Fault）
フーリエ変換・観測点データ： grflt12sx1-v3.f*** CHANGE OF SIGNS OF IMAGINARY PARTS OF FINAL RESULTS (FOR FFT) *Change sign (=1), or do not change sign (=0)
1*** OBSERVATION POINT ***NUMBER OF Obserbvation Points
1Z (m)
Observati Location: XY (m)163578-9182
0 : Luc→フーリエ変換の虚数部の符号を逆転（地震学定義→工学変換）
観測点：LUCの座標位置（X,Y,Z）① 表面波・位相速度の計算（分散曲線、Medium Response）• phs3sQ-v3.f
Medium Response of Love Wave (Up to Mode 5)

Dispersion Curves of Love Wave (Up to Mode 5)5000

4.50E-08
4.00E-084500

U 0 3.50E-08

4000

C 1

Medium Response

Phase and Group Velocities (m/s)

C 0

U 1

3500

C 2
3000

3.00E-08
M 02.50E-08

U 2

M 1

C 3 2.00E-08

M 2

U 3

2500

M 3

C 4 1.50E-08
U 4

2000

C 5
1500

位相速度

群速度

M 4
M 51.00E-08

U 5 5.00E-09

0.00E+00

1000
00.2

0.4

0.6

0.8

frequencies (Hz)Love波の分散曲線

100.2

0.4

0.6

0.8

1

frequency (Hz)Love波のMedium Response

H/Vスペクトル
1000

H/Vスペクトル（レイリー波）

100

H/V

10100.2

0.4

0.60.10.01

0.81M 0

M 1

M 2

M 3M 4M 5

M 6

Frequency (Hz)0.001

② 波数積分法による波形計算（振動数領域、速度波形出力）
・grflt12sx1-v２.f
修正表示定理（Hisada and Bielak, 2002）{}

U k (Y ; ω ) = ∫ Tik ( X , Y ; ω ) − TikS ( X , Y ) [Di ( X ; ω ) ]dΣΣ+ ∫ T ikS ( X , Y ) [D i ( X ; ω ) ]d ΣΣ動的グリーン関数（実体波・表面波）

静的グリーン関数（断層面近傍の項）

③ フーリエ逆変換（フィルター、時間領域）
・grfftspx.f
LUC（CEF断層のみ）
FPFNZ(UD)

速 度 (c m / s）

7550250-25010203040時間（秒）200

LUC（CEF断層のみ）

150

FPFN
Z(UD)

変 位 （ｃｍ ）100

100
5000-50-50

速度波形（出力波形）

1020

3040

時間（秒）

変位波形（速度波形を時間積分）

おわりに○強震動予測の必要性
標準地震波（告示波など）と様々な強震観測記録の違い（指向性パルス、フリングステップ）○波数積分法による強震動計算法
震源・伝播・サイト特性○計算コードと事例の紹介
計算コード・データのダウンロード
http://kouzou.cc.kogakuin.ac.jp/→ Open Dataへクリック