國立中央大學地震災害鏈風險評估及管理研究中心
大規模地震情境模擬對外資料

2021-東部地震情境模擬

部衝擊情境,預期成果係完成琉球隱沒帶情境震源參數擬定,並建立指標性特徵震源模型,依據擬定震源模型進行地震波模擬分析,最終提供災損分析使用之峰值地動加速度(peak ground acceleration, PGA)、峰值地動速度(peak ground velocity, PGV)、峰值地動位移(peak ground displacement, PGD)、0.3秒譜加速度(0.3-s of spectral acceleration, SA0.3)與 1.0 秒譜加速度(1.0-s of spectral acceleration, SA1.0),其地動數值涵蓋空間範圍包含整個臺灣東部主要地區(宜蘭、花蓮及台東北段)。
工作項目包含震源參數擬定、速度與構造模型建置、特徵震源模型建置、地震波模擬分析、不同數值模擬方法比較及地表地震動歷時與反應譜計算分析等工作。其中,第二章參數與方法章節,模擬分析所需參數蒐集及擬定章節,陳述如何針對琉球隱沒帶地震蒐集相關文獻用以擬定震源模型參數;地形與速度構造模型建置章節,說明供後續情境地動模擬使用且覆蓋分析區域之地震波速度構造模型;特徵震源模型建置章節,考慮琉球隱沒帶地震之特徵規模及臺灣東部地區可
能之地動影響模式,設定可能震源破裂模式,詳加擬定震源參數,並製表呈現;地震波模擬分析章節闡述使用方法;不同低頻震波數值模擬方法之比較。另,第三章則針對地表地震動歷時與反應譜計算分析,提出涵蓋研究區域內之模擬成
果。

第 2 章 參數與方法

        本章節闡述考量琉球隱沒帶情境地震而彙整之震源參數,並以此作為依據,以地震學、地體動力學與工程地震學角度,考慮對臺灣東部地區較具衝擊代表性之情境,擬定特徵震源模型。另對於速度、地形構造之數值模擬所需之模型建構進行考量。最終針對該單一特徵震源模型進行情境模擬地震動,模擬完整情境地震歷時,計算地震動參數,包含(1) 峰值地動加速度(PGA)、(2) 峰值地動速度(PGV)、(3) 峰值地動位移(PGD)、(4) 0.3 秒譜加速度(SA0.3)與(5) 1.0 秒譜加速度(SA1.0)等五項參數作為產出,以利後續災損推估小組推估臺灣東部地區震災情境。綜整上述說明,模擬分析程序詳見圖 2- 1。

圖 2- 1 震源情境地震動模擬分析流程
圖 2- 1 震源情境地震動模擬分析流程

2.1 震源參數擬定

        臺灣東北部外海琉球海溝具引致大規模隱沒帶地震之潛勢,歷史上於 1920年,此區域卽發生規模達 8.0 之地震事件,Theunissen et al. (2010)提出該區域潛勢斷層 ISZ:板塊介面孕震區域(interpolate seismogenic zone)模式為最可能發生之1920 年地震與琉球隱沒帶地震之破裂型態。此外,Hsu et al. (2012)使用 2005 年至 2010 年間之連續全球定位系统(Global Positioning System, GPS)資料,推估琉球隱沒帶斷層幾何及滑移虧損(slip-deficit),該區域具有發生震矩規模 7.5 至 8.7 的災害地震潛勢。本研究綜合考量上述論文研究成果,並根據國家地震工程中心SHACC level 3 計畫綜整近年臺灣海陸域斷層調查成果(PSHA SSHAC-3 SSC, 2019),作為擬定琉球隱沒帶情境震源參數之依據(圖 2- 2)。

圖 2- 2 SSHAC Level 3 琉球隱沒帶震源參數模型

2.2 地形與速度構造模型建置

        地震波傳模擬所需速度構造模型建置,係以臺灣地區速度構造研究(KuoChen et al., 2012)之三維速度分布數值,作為本研究模擬計算之模型建置基礎資料。該速度分布數據係藉由中央氣象局強震站及短週期地震站、中央氣象局及中央研究院寬頻地震站、大量布設臨時站所記錄到之大量地震波走時資訊,經由逆推 方法計算建構而成 。 參考三維速度分布數值所 強 調 之構造不均質(heterogeneity)與側向變化(lateral variation),提升了僅考量無方向性變化的一維速度模型可能影響地震模擬評估之不足。評估本先導研究情境模擬工作所需要之涵蓋空間範圍後,整合三維速度分布數值所建構之臺灣本島及外海區域速度構造模型如圖 2- 3 所示。

        此外,因數值方法特性,本研究所採行之地震動模擬方法於低頻模擬時可考慮地形起伏造成的震波散射效應,例如震波傳遞至山峰及山脊處時所造成振幅放大,山谷處造成振幅衰減(Lee et al., 2009)。因地形效應額外產生之表面波與散射波,預期會影響測站處記錄之波形。本研究計畫以解析度為一弧秒(約 1.85 km)之ETOPO1 (Amante and Eakins, 2009)做為初始地形模型(圖 2- 4),考量前述臺灣速度構造模型後,經穩定性測試後,計算出震波模擬適合之模型解析,將地形模型重採樣為 300 m 取一格點,作為數值計算網格水平方向格距,建立數值網格模型(圖 2- 5),模型空間範圍考量納入涵蓋琉球隱沒帶情境震源與臺灣東部目標場域之範圍。

圖 2- 3 臺灣本島及外海區域地下速度構造模型
圖 2- 3 臺灣本島及外海區域地下速度構造模型
圖 2- 4 臺灣本島及外海區域數值地形模型
圖 2- 4 臺灣本島及外海區域數值地形模型
圖 2- 5 臺灣本島及外海區域數值地形與速度構造整合模型
圖 2- 5 臺灣本島及外海區域數值地形與速度構造整合模型

2.3 特徵震源模型建置

        考量琉球隱沒帶具轉折之線形與延伸破裂面,本完整詳實地震動模擬須簡化其琉球隱沒帶線形,合理考量符合線型之斷層長度與寬度,進而假定特徵震源模型。

        特徵震源模型之參數擬定原則,則採用日本地震調查研究推進本部之 Recipe(Irikura and Miyake, 2011)方法,取其特徵震源擬定之標準化流程作為依據,針對琉球隱沒帶情境地震破裂,探討長週期地震動對臺灣東部區域之影響性,假定斷層面上具單一地栓(asperity)為主要能量釋放位置,將地栓中心置於沿斷層走向之中線,地栓上緣與斷層上緣重合,並考慮震源(即破裂起始點)位於地栓上緣最外海處,考量斷層由海床往深處破裂琉球隱沒帶情境地震震源模型以平面與側向示意,詳見圖 2- 6 與圖 2- 7;模型之震源位置(即破裂起始點)列於表 2- 1;斷層面之巨觀與微觀震源參數則分列於表 2- 2 與表 2- 3。

圖 2- 6 琉球隱沒帶情境震源之特徵震源模型平面示意圖
圖 2- 6 琉球隱沒帶情境震源之特徵震源模型平面示意圖
圖 2- 7 琉球隱沒帶情境震源之特徵震源模型之側向示意圖
圖 2- 7 琉球隱沒帶情境震源之特徵震源模型之側向示意圖
表 2- 1 琉球隱沒帶情震模型震源位置參數
表 2- 1 琉球隱沒帶情震模型震源位置參數
表 2- 2 琉球隱沒帶情震地震巨觀震源參數表
表 2- 2 琉球隱沒帶情震地震巨觀震源參數表
表 2- 3 琉球隱沒帶情震地震微觀震源參數表
表 2- 3 琉球隱沒帶情震地震微觀震源參數表

2.4 虛擬測站點位

2.4.1 點位分布設定

        本研究配合災損推估小組之評估網格設定,依據國家災害防救科技中心(National Science and Technology Center for Disaster Reduction, NCDR)格點解析度 500 公尺,確認分布範圍為東經 118.125 度至 122.005 度,北緯 21.895 度至 26.385 度,範圍涵蓋臺灣本島與離島範圍,共 132,712 網格點。圖 2- 8 為 NCDR 網格示意圖,黑色格點即為網格點位。

圖 2- 8 國家災害防救科技中心(NCDR)500 公尺網格點位
圖 2- 8 國家災害防救科技中心(NCDR)500 公尺網格點位
2.4.2 本研究設定點位

        對應災損推估小組後續評估工作,空間範圍與解析上,本研究基於 NCDR 500 公尺網格點擬定虛擬測站,進行地震動模擬計算。所採行點位分布範圍為東經120.500 度至 122.005 度,北緯 22.1 度至 25.5 度,共 24,755 點位,空間點位選取範圍系考慮琉球隱沒帶情境震源所在區域、臺灣東部生活圈,及關鍵基礎設施進而擬定之範圍。如圖 2- 9 所示,分布於臺灣東部陸地區域之灰色格點為 NCDR 500 公尺網格點位,視為虛擬測站;藍框區域為情境地震動模擬範圍。

圖 2- 9 本研究設定點位與震波模擬範圍對照圖
圖 2- 9 本研究設定點位與震波模擬範圍對照圖

2.5 地震波模擬分析

2.5.1 低頻模擬:三維有限差分法

        本研究之低頻地震波模擬方法,採用有限差分法(finite-difference method),於三維地形及地下速度構造模型設定中,透過擬定假想地震之情境(特徵震源模型)進行模擬,評估低頻地震波對特定空間點位、工址或測站產生之強地動影響。

        震波傳遞之模擬方法係基於物理機制下,考慮地震源與震波傳遞路徑中各種耦合效應,計算特定點位的震動歷時反應,其中震源項包含地震震源機制、斷層幾何與破裂型態,而路徑項則包含不同尺度的地下構造模型(包含波速、密度與衰減因子等)。以有限差分法進行震波模擬,主要係由求解質點運動方程與組構關係隨時間之變化,即為模擬質點受力後於空間中的運動行為與震波隨時空變化之傳遞過程。有限差分法為基於網格運算之數值方法,即為求網格系統中某格點位置之差分運算,運算過程僅需鄰近格點即可計算,因此可透過訊息傳輸介面(Message Passing Interface, Gropp et al., 1996)將數值方法改寫,並在叢集式電腦(PC cluster)進行平行化計算 (parallel computation),進而提高效率或進行更為高頻的模擬工作。

         本研究採用 Zhang et al. (2012)之牽引力鏡像有限差分法(traction-image finite difference method)可考慮(1) 三維地下構造與(2) 地表地形,進行震波模擬。該方法採用 MacCormack 格式(MacCormack, 1969; Gottlieb & Turkel, 1976)與 optimized Dispersion Relation Preserving (DRP/opt, e.g. Tam & Webb, 1993; Hixon, 1997)能有效降低數值頻散(numerical dispersion)效應。同時,該演算法採取完美匹配吸收層(perfectly matched layer, PML)做為吸收邊界條件(absorbing boundary condition, e.g.Berenger, 1994; Marcinkovich & Olsen, 2003),能有效降低人為設定模型邊界產生之反射波。

        為考慮地表地形或地下構造反射面效應,Zhang et al. (2012)透過坐標轉換,將地表地形與地下構造界面以曲線性坐標系統(curvilinear coordinate system)描述,進而使物理空間中各曲線性坐標位置透過 Jacobian 變換至電腦計算所使用的直角坐標系統 (圖 2- 10)。

圖 2- 10 座標轉換示例(摘自 Zhang et al., 2012)
圖 2- 10 座標轉換示例(摘自 Zhang et al., 2012)
2.5.2 高頻模擬:隨機式方法

        針對需求所擬定強地動模擬,須計算 PGA、SA0.3 等高頻地震動參數,模擬頻段須達 10 Hz,2.5.1 節之低頻模擬方法無法反應 1-10 Hz 頻段之高頻地震動時間歷時,故本研究計畫選用 EXSIM (Atkinson and Assatourians, 2015)隨機式方法進行高頻地震動模擬,補足人工合成地震歷時之高頻響應。

        EXSIM 係根據使用者給定之地震矩,以經驗公式給定斷層長寬形成一斷層面,而後離散化為一系列之子震源(subsource),其中的每一個子震源皆遵循點震源隨機模型(point-source stochastic model),依給定的地震歷時(duration)與剪力波輻射型態以隨機高斯噪訊(random Gaussian noise)模擬強地面加速度地動歷時。EXSIM 亦引入 Brune model,於傅氏譜(Fourier spectrum)考慮地震矩、應力參數、經驗衰減參數於頻率域之衰減效應,經快速傅立葉變換(fast Fourier transform,FFT)至時間域後,根據破裂速度疊加各子震源,並考慮場址效應修正函數,形成時間序列。本研究亦參考 Lee et al. (2015)根據場址 VS30 條件分類推估地震歷時之經驗公式,整合於高頻模擬計算。

         EXSIM 為相當具代表性之隨機式強地動模擬方法,除南加州地震中心(Southern California Earthquake Center, SCEC)已對 EXSIM 實行相關檢核(e.g. Atkinson and Assatourians, 2015)外,EXSIM 模擬結果亦已與隱沒帶地震事件檢核(e.g. Atkinson and Macias, 2009)。

        本研究以 EXSIM 模組進行隨機法強地動模擬時,給定琉球隱沒帶特徵震源模型之滑移量分布作為輸入震源模型,並針對模擬範圍,考慮適用臺灣區域之隨機式模型(Sokolov et al., 2006; Sokolov et al., 2009; Huang et al., 2017),詳細模擬參數列於表 2- 4。

表 2- 4 高頻模擬採用之隨機式模型參數表
2.5.3 混合法合成寬頻地震歷時

        考量前述 2.3 節地栓滑移量分布之隱沒帶震源情境特徵震源模型,並根據2.4.2 節所述採行之 24,755 網格點位,以三維有限差分法與隨機式分析方法分別模擬各點位之低頻與高頻地震波,並以剪力波波速輔以互相關(cross-correlation)方法修正各點位之震波到時(arrival-time)後疊加低頻與高頻地震波,以工程常用之 VS=760 m/s 基盤面做為地震動計算基準,以混合法合成寬頻地震歷時至該基盤面,後續便以該寬頻震動歷時進行 5 項地動參數計算,包含(1) 峰值地動加速度(PGA)、(2) 峰值地動速度(PGV)、(3) 峰值地動位移(PGD)、(4) 0.3 秒譜加速度(SA0.3)與(5) 1.0 秒譜加速度(SA1.0)。

        本研究以叢集式電腦進行平行運算,考慮單一情境,琉球隱沒帶情境地震之特徵震源震源模型共 9,016 子斷層與地表面 24,755 虛擬測站,其運算量共223,191,080 次點震源模擬,以 96 計算核心進行計算,低頻模擬約莫耗時 23 小時,高頻模耗時約 5.4 小時,混合法與場址修正耗時約 64 小時。後續若對單一斷層進行多次情境模擬,仍須審慎評估計算量與計算時間,合理使用並擴增計算資源使計算效益最大化。

2.5.4 不同數值模擬方法之比較

        許多種數值方法陸續發展應用於三維地震波模擬,例如上一章節(2.5.1)介紹的牽引力鏡像有限差分法(traction-image finite-difference method, Zhang et al.,2012),及其他常見方法如,交錯式網格有限差分法(staggered-grid finite-difference method, Open-source Seismic Wave Propagation Code, OpenSWPC, Maeda et al.,2017) 和 譜 元 素 法 (spectral-element method, Komatitsch and Vilotte, 1998 ;Komatitsch and Tromp, 1999)等方法。

        關於不同數值方法之間的量化比較工作,第三屆的地表地質效應研討會(International symposium on the effects of surface geology, ESG 2006)已有團隊(Chaljub et al., 2010)進行詳盡的測試。該團隊比較四種不同數值方法,包括任意高階導數非連續伽遼金法(arbitrary high-order derivative discontinuous Galerkin
method, ADER-DGM, Käser et al., 2006) ,有限差分法(finite-difference method, FDM, Kristek et al., 2009),譜元素法(spectral-element method, SEM, Chaljub, 2009),以及第二種譜元素法(spectral-element method, SEM, Stupazzini, 2009 and Stupazzini et al., 2009)。結論認為,雖然現今數值方法仍在發展階段,但是在適當地使用下,四種方法皆可以正確處理複雜的三維速度模型異質性,獲得可靠的地動預估結果。此外,也建議採用兩種以上相對準確的方法進行模擬,例如有限差分法以及譜元素法,如此可以獲得更加準確可靠的結果。

        本研究使用牽引力鏡像有限差分法(以下簡稱為 FDM1)、交錯式網格有限差分法(以下簡稱為 FDM2)以及譜元素法(以下簡稱為 SEM),進而比較三種數值方法採用相同參數時所模擬之合成波形,透過量化合成波形差異,理解不同的數值方法及不同的操作過程所可能出現的差異。比較過程中共測試兩組模型,第一組模型為速度均勻不具側向變化,且地形平坦之簡化模型,模型剪力波速為 2.0 km/s,壓縮波速為 3.5 km/s,密度為 2.34 g/cm3;採用點震源,震源走向 90°,傾角 90°,滑移角 180°,震源深度 8.0 公里位於模型原點(0,0);模型尺度範圍及測站配置如圖 2- 11 所示。第一組模型採用三種數值方法進行比較,三種方法之合成波形比較,以 LINE_001 測站陣列為例,如圖 2- 12 所示。 根據震源機制,理論上 X 與Z 分量振幅為 0,僅 Y 分量有值存在;三種方法之合成波皆呈現非常相近的結果,包含振幅與波形相位。

         而為了量化合成波形之間的差異,本研究導入 Kristekova et al. (2009)所提出的時間頻率域擬合度 (Time-frequency misfit) 評估方法,以 TF_MISFIT_GOF_CRITERIA 套件,對合成波形之間進行差異的量化計算。Kristekova et al. (2009)對波形之間的差異在頻率域與時間域中進行量化,並給出量化分數介於 0~10 之間,當擬合度高於 8.5 時屬於”excellent”的擬合度,擬合度高於 6.5 時屬於”good”的擬合度,擬合度高於 4.5 時屬於”fair”的擬合度,擬合度低於 4.5 時屬於”poor”的擬合度。本研究皆以 FDM1 之合成波形做為參考波形,分別對 FDM2及 SEM 合成波形進行時間頻率域擬合度計算。圖 2- 13 為 FDM2 與 FDM1 之相位擬合度及振幅擬合度,其中相位擬合度只有在測站接近震源時,擬合度略為下降,或位於震源上方(Station ID=41)時之奇異點大幅下降,而其他測站之波形相位擬合度皆非常高;而振幅擬合度則呈現 Y 分量整體略為下降,然而平均擬合度仍在 8.5 分以上,X 與 Z 分量則皆趨近滿分。圖 2- 14 為 SEM 與 FDM1 之相位擬合度及振幅擬合度,其中相位擬合度只有在震源上方(Station ID=41)時之奇異點大幅下降,其餘測站皆趨近滿分擬合度;而振幅擬合度也呈現 Y 分量整體略為下降,平均擬合度高於 9 分,X 與 Z 分量則皆趨近滿分。

        根據測試過程得知,FDM2 與 FDM1 之合成波形差異來源主要包含兩點:第一點為非交錯式網格有限差分法與交錯式網格有限差分法之數值方法本身對於座標系統之定義所帶來的誤差,即 FDM2 之實際測站位置與 FDM1 差距為半個網格點,以本測試來說大約差 12.5 公尺;第二點為震源時間函數定義不同,FDM1 採用 GAUSSIAN 函數,而 FDM2 採用 HEERMAN 函數做為震源時間函數。而 SEM 與FDM1 之合成波形差異比例極低,此差異來源屬不同的數值方法間的正常計算誤差;惟仍須注意 SEM 與 FDM 模型網格化時可能產生之座標誤差。

        測試之第二組模型為考慮側向速度變化及地形起伏之模型,速度模型採用(Kuo-Chen et al., 2012)之三維速度分布數值,而地形起伏採用 ETOPO1 (Amante and Eakins, 2009)做為初始地形模型。震源則採用前期報告之山腳斷層南段破裂有限斷層模型。在第二組測試中僅以 FDM1 與 SEM 進行比較。圖 2- 15 為選取的 6 個測站與斷層的相對位置,而經選取 6 個測站之 FDM1 與 SEM 合成波形比較如圖2- 16 所示。合成波形比較顯示,在多數測站中,FDM1 與 SEM 之合成波形皆十分相似(圖 2- 16a、b、c、d);僅在近震源處(圖 2- 16e、f)有極少數測站在振幅上出現較明顯的差異。本研究在空間上均勻選取 1492 個測站進行時間頻率域擬合度計算,以量化評估 FDM1 與 SEM 之合成波形差異,圖 2- 17 為振幅擬合度之分布,整體平均之擬合度達”excellent”,僅少數測點屬”good”。圖 2- 18 為相位擬合度之分布,整體平均之擬合度均達”excellent”,波形相位擬和極佳。本研究研判,該較大波形差異測點的產生,可能與本研究半盲測試過程中,模型(速度,地形,座標)仍有些微不同所致。

        經由兩組模型測試指出,有限差分法與譜元素法在適當的操作下,其計算結果非常接近,唯須注意操作過程的各項參數之設計,可能影響計算結果產生差異。

圖 2- 11 第一組簡化模型之尺度與震源及測站之配置
圖 2- 11 第一組簡化模型之尺度與震源及測站之配置
圖 2- 12 第一組簡化模型三種數值方法之合成波形比較,(a)(b)(c)為牽引力鏡像有限差分法之三分量合成波形;(d)(e)(f)為交錯網格有限差分法之三分量合成波 形;(g)(h)(i)為譜元法之三分量合成波形;所有波形振幅以牽引力鏡像有限差分 法之 Y 分量最大振幅進行正規化;波形經低通濾波 2Hz
圖 2- 12 第一組簡化模型三種數值方法之合成波形比較,(a)(b)(c)為牽引力鏡像有限差分法之三分量合成波形;(d)(e)(f)為交錯網格有限差分法之三分量合成波 形;(g)(h)(i)為譜元法之三分量合成波形;所有波形振幅以牽引力鏡像有限差分 法之 Y 分量最大振幅進行正規化;波形經低通濾波 2Hz
圖 2- 13 FDM2 與 FDM1 之相位擬合度及振幅擬合度
圖 2- 13 FDM2 與 FDM1 之相位擬合度及振幅擬合度
圖 2- 14 SEM 與 FDM1 之相位擬合度及振幅擬合度
圖 2- 14 SEM 與 FDM1 之相位擬合度及振幅擬合度
圖 2- 15 挑選之 6 個測站(紅色點)及斷層面位置分布(深藍與暗紅),黑色點為模擬使用之所有測站位置,共 13216 個站
圖 2- 15 挑選之 6 個測站(紅色點)及斷層面位置分布(深藍與暗紅),黑色點為模擬使用之所有測站位置,共 13216 個站
圖 2- 16 挑選 6 個測站之 FDM1 與 SEM 合成波形比較
圖 2- 16 挑選 6 個測站之 FDM1 與 SEM 合成波形比較
圖 2- 17 SEM 與 FDM1 合成波形三分量之振幅擬合度分布
圖 2- 17 SEM 與 FDM1 合成波形三分量之振幅擬合度分布
圖 2- 18 SEM 與 FDM1 合成波形三分量之相位擬合度分布

第 3 章 地震動歷時與反應譜計算分析

        依前述地震波模擬分析方法,考慮琉球隱沒帶情境地震之特徵震源模型,分別進行低頻與高頻模擬,進而以混合法合成於工程基盤之寬頻加速度地震歷時,爾後分別積分一次得速度歷時、積分二次得位移歷時,各地震歷時之最大振幅處即分別為(1) 峰值地動加速度(PGA)、(2) 峰值地動速度(PGV)、(3) 峰值地動位移(PGD)等地動數值。另一方面,選取加速度地震歷時作為輸入地震歷時,考慮 5%阻尼計算該地震歷時之反應譜(response spectrum),爾後分別取週期 0.3 秒與週期 1.0 秒之譜加速度作為反應譜地動數值。本研究所模擬之地震動數值敘述如下:

        考慮震源位於地栓上緣最外海處,斷層由海床往深處破裂,因斷層破裂的方向指向臺灣本島花蓮方向,破裂的方向性效應強烈,此起情境引致臺灣花蓮區域地震動甚鉅,PGA 分布(圖 3-1a)顯示最大值為 594.4 cm/s2,出現在花蓮縣豐濱鄉一帶;PGV 分布如圖 3-1b,最大值為 141.5 cm/s,僅在非常少數區域出現此極端值,而豐濱鄉、瑞穗鄉、光復鄉及鳳林鎮區域皆有超過 80 cm/s 的 PGV 分布;該區域附近具有最大對應之震度(圖 3-2a)為 6 強;PGD 分布(圖 3-2b)在秀林鄉、萬榮鄉南部及豐濱鄉沿海一帶有超過 75cm 的分布,而最大值達到 98cm;SA0.3 分布(圖 3-3a)與 SA1.0 分布(圖 3-3b)之最大值則分別為 1.01 g 與 0.73 g。

圖 3-1 琉球隱沒帶情境地震模擬臺灣東部地區 PGA 分布(圖 a)與震度分布(圖 b)
圖 3-1 琉球隱沒帶情境地震模擬臺灣東部地區 PGA 分布(圖 a)與震度分布(圖 b)
圖 3-2 琉球隱沒帶情境地震模擬臺灣東部地區 PGV 分布(圖 a)與 PGD 分布(圖 b)
圖 3-2 琉球隱沒帶情境地震模擬臺灣東部地區 PGV 分布(圖 a)與 PGD 分布(圖 b)
圖 3-3 琉球隱沒帶情境地震模擬臺灣東部地區之 SA0.3 分布(圖 a)與 SA1.0 分布(圖 b)
圖 3-3 琉球隱沒帶情境地震模擬臺灣東部地區之 SA0.3 分布(圖 a)與 SA1.0 分布(圖 b)

第 4 章 結論與建議

        過去,琉球隱沒帶地震事件曾引致地震災害,如 1920 年琉球隱沒帶地震事件、1986 花蓮外海地震事件等,皆造成傷亡與建物毀損。本研究計畫擇定琉球隱沒帶作為假想目標,考慮單一地栓位置與震源向西破裂起始位置、震波傳遞等因素,擬定單一震源破裂情境,針對臺灣東部地區以混合法進行地震動模擬,完成情境模擬地動參數,包含峰值地動加速度、峰值地動速度、峰值地動位移、0.3 與1.0 秒譜加速度等,供後續災損推估小組進行評估。研究成果顯示,該震源模式對臺灣東部區域引致高地動量值,尤其是花蓮縣豐濱鄉、瑞穗鄉、光復鄉及鳳林鎮附近區域,具有 6 強以上震度;且引起之地動值兼具短周期及長周期之地震動,可能引致較大規模之地震災害。

        然而,本研究僅考慮單一地栓震源破裂情境,倘若琉球隱沒帶破裂,不全然為本研究所假想之震源破裂型態以及推估規模。故本研究報告建議如下:

  1. 針對單一地栓震源破裂面,未來可考慮不同地栓位置之滑移量分布模式、震源破裂位置等多種破裂情境,對各情境產出之地震動參數進行統計分析,有效量化地震動中值與極值,估算地震動範圍。
  2. 對於分析地點進行場址放大因子研究,提供更高空間解析之準確場址放大,藉以提升地動空間解析。
  3. 有限差分法或譜元素法皆為低頻三維地動模擬有效工具,透過量化比較指出,兩者計算結果相似;但須注意操作過程使用之參數及流程之正確性。

參考文獻

Amante, C. and Eakins, B. W. (2009). ETOPO1 1 arc-minute global relief model:
         procedures, data sources and analysis, NOAA Technical Memorandum, NESDIS NGDC-24, 19 pp.
Atkinson, G. M and Assatourians, K. (2015). Implementation and validation of EXSIM
         (a stochastic finite-fault ground-motion simulation algorithm) on the SCEC broadband platform, Seismol. Res.. Lett., 86(1), 48–60, DOI:10.1785/0220140097.
Atkinson, G. M. and Macias, M. (2009). Predicted ground motions for great interface
         earthquakes in the Cascadia subduction zone, Bull. Seism. Soc. Am., 99(3), 1552–1578.
Berenger, J. A., (1994). Perfectly matched layer for the absorption of electromagnetic
         waves, J. Comput. Geophys., 114, 185–200.
Chaljub, E. (2009). Spectral element modeling of 3D wave propagation in the Alpine
          valley of Grenoble, France, in ESG 2006, Third Intl. Symposium on the Effects of Surface Geology on Seismic Motion, P.-Y. Bard, E. Chaljub, C. Cornou, F. Cotton, and P. Guéguen (Editors), LCPC Editions, ISSN 1628-4704, Vol. 2, 1467–1473.
Chaljub, E., Moczo, P., Tsuno, S., Bard, P. Y., Kristek, J., Käser, M., … & Kristekova, M.
          (2010). Quantitative comparison of four numerical predictions of 3D ground motion in the Grenoble Valley, France. Bulletin of the Seismological Society of America, 100(4), 1427-1455.
Gottlieb, D. and Turkel, E. (1976). Dissipative 2-4 methods for time-dependent
           problems, Math. Comput., 30, 703–723.
Gropp, W., Lusk, E., Doss, N., and Skjellum, A. (1996). A high-performance, portable
           implementation of the MPI message passing interface standard, Parallel Comput., 22(6), 789–828.
Hixon, R. (1997). On increasing the accuracy of MacCormack schemes for aero acoustic
           applications, AIAA Paper, doi: 10.2514/6.1997–1586.
Hsu, Y. J., Ando, M., Yu, S. B., and Simons, M. (2012). The potential for a great
            earthquake along the southernmost Ryukyu subduction zone, Geophys. Res. Lett.,39, L14302, doi:10.1029/2012GL052764.
Huang J. Y., Wen, K. L., Lin, C. M., Kuo, C. H., Chen, C. T., and Chang, S. C. (2017). Site
            correction of a high-frequency strong-ground-motion simulation based on an empirical transfer function, J. Asian Earth Sci., 138, 399–415.
Irikura, K. and Miyake, H. (2011). Recipe for predicting strong ground motion from
             crustal earthquake scenarios, Pure Appl. Geophys., 168, 85–104,doi:10.1007/s00024-010-0150-9.
Käser, M., M. Dumbser, and J. de la Puente Alvarez (2006). An efficient ADER-DG
             method for 3-dimensional seismic wave propagation in media with complex geometry, in ESG 2006, Third Intl. Symposium on the Effects of Surface Geology on Seismic Motion, P.-Y. Bard, E. Chaljub, C. Cornou, F. Cotton, and P. Guéguen(Editors), LCPC EditionsVol. 1, pp. 455–464, Grenoble.
Komatitsch, D., & Vilotte, J. P. (1998). The spectral element method: an efficient tool
             to simulate the seismic response of 2D and 3D geological structures. Bulletin of the seismological society of America, 88(2), 368-392.
Komatitsch, D., & Tromp, J. (1999). Introduction to the spectral element method for
             three-dimensional seismic wave propagation. Geophysical journal international,
139(3), 806-822.
Kristeková, M., Kristek, J., & Moczo, P. (2009). Time-frequency misfit and goodness-of fit criteria for    quantitative comparison of time signals. Geophysical Journal International, 178(2), 813-825.
Kristek, J., P. Moczo, and P. Pazak (2009). Numerical modeling of earthquake motion in
            Grenoble basin, France, using a 4th-order velocity- stress arbitrary discontinuous staggered-grid FD scheme, in ESG 2006, Third Intl. Symposium on the Effects of Surface Geology on Seismic Motion, P.-Y. Bard, E. Chaljub, C. Cornou, F. Cotton,
and P. Guéguen (Editors), LCPC Editions, ISSN 1628-4704, Vol. 2, p. 1517–1526.
Kuo-Chen, H., Wu, F. T., and Roecher, S. W. (2012). Three-dimensional P velocity
            structures of the lithosphere beneath Taiwan from the analysis of Taiger and related seismic data sets, J. Geophys. Res., B6.
Lee, S. J., Komatitsch, D., Huang, B. S., and Tromp, J. (2009). Effects of topography on
             seismic-wave propagation: an example from northern Taiwan, Bull. Seism. Soc. Am., 99, 314–325.
Lee, Y. T., Ma, K. F., Wang, Y. J., and Wen, K. L. (2015). An empirical equation of effective
             shaking duration for moderate to large earthquakes, Nat. Hazards, 75(2), 1779–1793.
MacCormack, R. W., (1969). The effect of viscosity in hypervelocity impact cratering,
              AIAA Paper, 69–354.
Maeda, T., Takemura, S., & Furumura, T. (2017). OpenSWPC: an open-source integrated
               parallel simulation code for modeling seismic wave propagation in 3D heterogeneous viscoelastic media. Earth, Planets and Space, 69(1), 1-20.
Marcinkovich, C. and Olsen, K. (2003). On the implementation of perfectly matched
                layers in a three-dimensional fourth-order velocity-stress finite difference scheme, J. Geophy. Res., 108, doi: 10.1029/2002JB002235. issn: 0148-0227.
PSHA SSHAC-3 SSC, 2019, Development of the Hazard Input Document for Taiwan
               using SSHAC Level 3 Methodology. Volume 2: SSC Technical Report, Taiwan Power Company.
Sokolov, V., Loh, C. H., and Jean, W. Y. (2006). Strong ground motion source scaling and
               attenuation models for earthquakes located in different source zones in Taiwan,4th International Conference on Earthquake Engineering, , Paper No. 003.
Sokolov, V., Wen, K. L., Miksat, J., Wenzel, F., and Chen, C. T. (2009). Analysis of Taipei
                basin response for earthquakes of various depths and locations using empirical data, Terr. Atmos. Ocean. Sci., 20, 687–702, doi: 10.3319/TAO.2008.10.15.01(T).
Stupazzini, M. (2009). 3D Ground Motion Simulation of the Grenoble Valley by
                GeoELSE, in ESG 2006, Third Intl. Symposium on the Effects of Surface Geology on Seismic Motion, P.-Y. Bard, E. Chaljub, C. Cornou, F. Cotton, and P. Guéguen(Editors), LCPC Editions, ISSN 1628-4704, Vol. 2, 1551–1560.
Stupazzini, M., R. Paolucci, and H. Igel (2009). Near-fault earthquake ground motion
                simulation in the Grenoble Valley by a highperformance spectral element code,Bull. Seismol. Soc. Am. 99, 286–301.
Tam, C. and Webb, J. C. (1993). Dispersion-Relation-Preserving finite difference
               schemes for computational acoustics, J. Comput. Phys., 107, 262–281.
Theunissen, T., Font, Y., Lallemand, S., and Liang, W. T. (2010). The largest
               instrumentally recorded earthquake in Taiwan: Revised location and magnitude, and tectonic significance of the 1920 event, Geophys. J. Int., 183, 1119–1133, doi:10.1111/j.1365-246X.2010.04813.x.
Zhang, W., Zhang, Z. G., and Chen, X. F. (2012). Three-dimensional elastic wave
               numerical modelling in the presence of surface topography by a collocated-grid finite-difference method on curvilinear grids, Geophys. J. Int., 190, 358–378.
行政院 (2016) 民國 105 年災害防救白皮書。
中央研究院 (2015) 大規模地震災害防治策略建議書。中央研究院報告 No.13