物理地震危險性分析中的地震動模擬方法研究進展

李 奇，陳學良，高孟潭，李宗超

（中國地震局地球物理研究所，北京 100081）

引言

地震災害是全世界都要面臨的嚴峻問題，嚴重威脅人民的生命財產安全。破壞性地震在全世界范圍內都曾發生過，如1906 年的舊金山大地震、1923 年的關東地震、1976 年的唐山地震、1999 年的臺灣集集地震、2008年的汶川地震和2016年的熊本地震［1－9］等，這些地震給全世界造成了重大生命財產損失。盡管工程建筑水平以及抗震措施相對于以前已經有了很大的改善，但是破壞性地震依然會給國家帶來嚴重的損失。

因此，防震減災工作任重道遠。造成地震災害的原因除了抗震措施不足外，主要是人們缺乏對地震的清楚認知，未能夠正確設防，即使在地震研究充分的地區，仍然沒能很好地解決。雖然地震難以精確預測，但能通過歷史地震以及地質構造特點來推斷區域未來可能發生的地震動特點和最大震級，這種方法稱為地震危險性分析。地震危險性分析是評定某一區域內工程在有效使用期內遭受一定地震動強度的危險性，分析結果并以設防烈度和設計地震動表示，為工程抗震設計、震害預防和城市規劃提供地震動輸入［10－11］。由此可見：地震危險性分析是城市規劃和抗震設防的重要手段，能夠為抗震御災提供科學的依據。在面對城市遭受過破壞性地震的慘痛背景下，我國極為重視并在21世紀初全面開始和研究大城市的活動斷裂探測和地震危險性評價，為我國的經濟發展和社會穩定打下了堅實的基礎。

由于地震波傳播過程中存在非常復雜的隨機性和不確定性，但又必須了解和研究地震動的特點，因此地震工程學家需要深入研究地震動特點來為建筑物等結構提供結構抗震設計分析［12－18］。通過地震觀測數據來研究地震動特點是常用方法，但在實際分析應用上仍然存在一定的局限性，而隨著數值模擬發展越發成熟以及計算條件越發優越，數值模擬已成為研究地震動重要方法之一。地震動數值模擬應用廣泛，可用于抗震設計、強地震動基本參數標定、建筑物震害分析、地震動參數輸入以及地震危險性分析等［19－22］，本文主要敘述地震動數值模擬在地震危險性分析中的應用。

目前，工程地震大多使用的是經驗校準的地震動預測方程來量化未來地震的震動強度，并將得到的地震動作為工程的地震動輸入［23－24］。經驗方法雖然成熟，但是在沒有觀測記錄的情況下，經驗方法是有局限性的，如深沉積盆地的影響［25－26］，而地震動數值模擬是有可以突破這些限制。SCEC/CME就利用地震動數值模擬來預測未來可能發生的地震動［27－30］。這種基于物理實際和真實速度結構的地震動數值模擬［31－32］，能夠很好的反應地層速度結構、斷層的幾何真實異質性和盆地效應等影響［25，33－36］，雖然處在發展完善階段，但其已經發揮了重要作用，在沒有很好觀測條件的情況下，能夠很好地校正經驗方法并得到相對更加合理的結果。不過這種方法同樣也面臨著重大挑戰，如最大截至頻率為1Hz、如何將研究區域擴大和如何設置凹凸提參數等。南加州地震中心（SCEC）的CyberShake 項目利用真實斷層構建的震源模型和地區三維速度結構模擬了未來地震動并繪制地震危險性曲線，相對于傳統地震動預測方程，CyberShake能夠更好地模擬出方向性效應和反應盆地效應等波的傳播特征，能夠為特定場點提供更加可靠的設防依據。

隨著計算科學和地震科學的進一步發展，大量的地震數據處理和地震動模擬不再像以前那樣困難。除了美國的CyberShake 項目，還有日本的Recipe 方法等，都在嘗試用基于物理的地震動數值模擬的方法來計算地震危險性。目前來看：基于地震動數值模擬的地震危險性分析雖然存在不足，但是不能忽略其在未來的發展潛力。本文綜述了前人的研究成果，舉例介紹了正在進行的研究，探討未來該方法的研究趨勢。

1 地震危險性評價方法

地震危險性問題是對研究區域的在未來可能發生的地震引發的地震動評估，這種地震動危險性分析目前最常用的有兩種方法：一是概率性方法是以概率的方式分析、評價與表述未來地震危險性；二是確定性方法由給定的距離和震級下的設定地震引起的地震強度。確定性方法包括歷史地震法和地震構造法，取結果之大者作為地震動危險性確定性分析的結果。

1.1 概率性地震危險性分析方法

早在20 世紀60 年代，CORNELL［23］就提出了概率性地震危險性評價方法（PSHA）和潛在震源的概念，并且分點源、線源和面源探討了其他地震危險性的計算方法。自從CORNELL 提出PSHA 以來，就在世界范圍內得到了廣泛應用，包括地震區劃圖的繪制、抗震設計規范的應用和核電設施的安全性評價等。20 世紀70年代，美國繪制了基于PSHA的全國性地震動區劃圖。隨著概率性地震危險性不斷發展和完善，不同學者也做了相關研究。KIUREGHIAN 等［37］提出了大地震的破裂長度對計算的危險性有顯著影響；REITER［38］把PSHA 大致分為四個環節：（1）劃分潛在震源區；（2）確定潛在震源區地震活動性參數；（3）估計地震影響；（4）計算場地的地震危險性，通過這四個環節來計算地震危險性。

PSHA 引入中國后，高孟潭［39］對PSHA 進一步研究并改進已讓其符合中國國情—CPSHA，基于該方法編制了“中國地震烈度區劃圖（1990）”。但是概率性地震危險性分析方法存在諸多的不確定性等問題［40－49］，如何更好地完善PSHA是學者們不斷研究和探究的重點。

1.2 確定性地震危險性方法

確定性地震危險性方法（DSHA）分為兩類：歷史地震法和地震構造法。歷史地震法是基于歷史資料和調查資料選取最大地震動。地震構造法是基于區域內大地震的構造條件，評估區域內最大潛在地震。由于DSHA 是基于地震科學基礎，因此在實際地震危險性分析中也得到了廣泛應用。美國在編制新一代的地震區劃圖時用確定分析的結果修正PSHA 結果。因此基于DSHA 的危險分析方法在抗震御災上發揮著重要作用。WANG等［50］用DSHA繪制臺灣地區的地震動危害圖，發現50年超越概率10%的平均PGA（峰值加速度）和PSHA 計算結果基本一致；ZUCCOLO 等［51］將DSHA 和PSHA 預測結果和真實地震活動做了對比，發現DSHA 預測結果更加貼近真實記錄；MUALCHIN［52］和MURMAN 等［53］建議在地震設計中使用DSHA。所以即使在PSHA占據主要評定方法的情況下，DSHA的作用也不能忽視。

2 地震動數值模擬

豐富的地震記錄能夠在研究地震時起到關鍵作用，但是對于地震記錄稍顯不足或稀缺的地區（包括地震活動性處在相對較弱時期的地震帶），越發成熟的數值模擬是彌補這些區域地震危險性研究的必要方法。目前地震動數值模擬方法主要分為：確定性方法、隨機性方法和混合方法。

2.1 確定性方法

確定性方法是在基于AKI等［54］的彈性動力學表示定理，將地震動表達為震源函數與格林函數的時空卷積。目前確定性方法主要包括，有限元方法、有限差分法、離散波數法和譜元法等［55－66］。眾多的學者用確定性方法研究了斷層的盆地效應和滑沖效應、速度大脈沖等地震動特征［67－73］。

對于長周期地震動（＜1Hz），確定性方法能夠很好的模擬，然而對于短周期分量，確定性方法表現較差，這是因為在短周期地震動包含較強的隨機性，地殼模型和破裂過程中存在不均勻性，因此單一固定的震源時間函數不能反映隨機性和不均勻性。

2.2 隨機性方法

隨機性方法是基于高斯帶限白噪聲的隨機振動理論方法［74－75］，或者是基于小震的經驗格林函數方法［76－78］，該方法適用于短周期（頻率一般大于1Hz）地震動。前者在小震和遠距離觀測的情況下考慮為點源模型［79］，在大震和近距離時考慮為有限斷層模型［80－81］。經驗格林函數法是將研究區域符合要求的（主震震源機制和信噪比高的小震記錄等）地震記錄作為經驗格林函數模擬地震動［82－83］，經驗格林函數法已經包含了場地效應和路徑效應等信息，因此不用再計算復雜的理論格林函數。經驗格林函數的不足之處是難以模擬小震記錄匱乏的地區。

2.3 混合方法

混合方法是將確定性方法和隨機性方法結合起來利用各自的優勢，在高頻段采用隨機性方法，在低頻段用確定性方法，在分別模擬后經過適當的高通和低通濾波之后，在時域中疊加合成寬頻帶的強地震動［84－88］。雖然混合方法利用了確定性方法和隨機性方法的優點，但是混合方法也存在一定的局限性［89］，比如確定性信號結束與隨機信號開始的銜接處容易產生反應譜的畸變，而這段頻率往往是工程上關注的頻段范圍，因此混合方法還需要不斷地優化。

3 地震動數值模擬在物理地震危險性分析的應用

3.1 CyberShake項目

OWENS［90］提出將計算科學和地震科學結合起來之后，南加州地震中心（SCEC）的CME（Community Modeling Environment）項目頗受關注［29］。CME 項目旨在系統地研究地震科學，其重要目標就是研究區域地震危險性并為未來工程建設和震害防御提供科學依據。初期項目開展后，SCEC 分為兩個研究團隊。前者開發基于Java程序的算法—OpenSHA［91］。OpenSHA 能實現多種地震危險性概念計算包括ERF（地震預測模型）、IMR（烈度度量關系）和IMT（烈度度量類型），OpenSHA 能為CME 提供了更加靈活的建模環境，讓其在復雜計算需求下仍能適用。后者構建南加州的三維速度模型（CVM-V3.0）［92－93］，其能反應真實地層的起伏并為地震動數值模擬提供更高的計算精度。模擬地震動與觀測記錄的對比表明精細的三維速度結構能更好地模擬地震動的傳播。

烈度度量關系（IMRs）和地震破裂預測模型（ERFs）是計算PSHA 的兩個關鍵。ERF 是在給定確定場點條件下，能夠生成該區域可能發生的地震集合。此外，還能給出每一個地震的震級和發生概率等。這些地震的集合是根據歷史地震、斷層活動性和斷層尺度等來確定的。IMR 是地震波隨距離的衰減關系。起初，IMR 使用的是經驗衰減關系，后來基于波形的IMR 發展起來（即后來的CyberShake 項目），SCEC 科學家們希望CyberShake能夠取代傳統的衰減關系來推測未來南加州地區的地震危險性。

近年來，SCEC/CME 還開展了多項研究。面向工程的地震動模擬寬頻帶平臺（Broadband Plantform）［94－98］，基于物理的全概率地震危險性分析項目——CyberShake［99－100］，高頻帶擴展項目（High-F）［101］，基于三維層析成像技術的三維速度結構模型（F3DT）等［31，102］。其中：CyberShake將確定性震源和波傳播效應納入地震動模擬中，由此來計算地震危險性，計算流程見表1。為了更好地構建南加州地震實驗場的地震危險概率分析的精度和提高模擬的頻率上限，CyberShake 也在不斷的發展和完善，如ERF 的更新換代和構建更加精細的CVM等。

表1 CyberShake的計算流程［99］Table 1 The workflows of CyberShake

CyberShake 計算地震危險性時，采用的是確定性震源，大量的計算不可避免，計算范圍內的地震震級越大，產生的破裂變化也越多，因此計算會消耗大量的計算內存和時間。為了能夠降低計算成本，利用互易定理（reciprocity）產生大量的格林應變張量（SGT）來合成地震圖從而節約計算成本［34，103－104］。雖然利用互易定理降低了計算成本，但是互易定理的計算要求是彈性介質［105－107］，難以反應非彈性非線性的性質。為了能將非彈性的性質考慮進地震動模擬中，將近斷層的非彈性性質納入到破裂模型的震源效應里，將近地表的非彈性性質納入到路徑效應里，這樣即能保證其彈性介質要求又能考慮非彈性影響。SCEC 還采用了一種新的范式（paradigm）來取代傳統震源效應、路徑效應和場地效應的解耦計算［108－109］，在新的范式下，地震動模擬可以計算非線性響應，并且震源、路徑和場地效應是耦合計算的，通過新的paradigm 能夠更加恰當地體現斷層外的非彈性性質和淺表地層的流變性質。

不過，CyberShake 項目仍然具有局限性：①波動場的計算模型是偽動力學（運動學）震源模型，并非真正意義上的動力學模型；②相對于工程感興趣的頻段（1～10 Hz），項目計算的最大截止頻率為1 Hz；③彈性計算條件限制，無法實現真正意義上的非彈性計算條件。這些都是需要面對的巨大難題，相信隨著認識的增加和計算效率的提升，CyberShake會不斷地完善和發展有望在未來地震危險性分析中扮演重要角色［110－115］。

確定性震源的地震動模擬方法能夠模擬出盆地效應和速度大脈沖等特征，能反應特定地質構造下的地震動特征，這是衰減關系難以做得到，但是CyberShake 目前還處在發展階段存在限制和不足如：地層速度結構不夠精細和斷層系統的理解局限等，未來會隨著計算能力的加強和認識探索的不斷深入，CyberShake極具潛力。

3.2 Recipe方法

ERC（Earthquake Research Committee）在2005 年出版了兩種日本國家地震災害圖［116］：一種是基于PSHA的地震災害圖；另一種則是指定震源斷層的地震災害圖。前者是表明全國不同地區可能發生的強烈地震的可能性，而后者想要預測區域內各個地震共同作用下的強地震動，ERC 把這種預測強地震動的方法稱為Recipe 方法。Recipe 方法是由IRIKURA［117－118］提出的，其將地震動的長周期和短周期分量分別用有限差分法和隨機格林函數法計算，然后借此來預測由活斷層引發的強地面運動。Recipe 方法的關鍵點在于震源模型的建立，IRILURA 把震源參數分為三類：外部震源參數（outer fault parameters）、內部震源參數（inner fault parameters）和額外震源參數（extra fault parameters）。

外部斷層參數包含整個斷層的幾何尺寸和地震矩。參數可由下式計算［119－120］：

式中：M0為地震矩；為平均應力降；S為斷層面積；W和L分為是斷層的寬度和長度；a和b分別取1.4 ×10-2和1。式（1）適用于地震矩小于1019Nm的情況和式（2）用于大于1019Nm的情況。

內部斷層參數有凹凸體參數（包括面積、應力降等）和加速度震源譜幅值（amplitude of acceleration source spectrum level）。內部斷層參數表征斷層非均質性的特點。參數可由下式計算［121－123］：

式中：Δσa為凹凸體上的應力降；Sa為凹凸體面積；A0和Aa0分別是加速度震源譜幅值和凹凸體上的加速度震源譜幅值；β和vr分別是介質中的S波傳播波速和破裂速度。對于內陸地震SaS取0.22，對于俯沖型地震SaS取0.25。

額外斷層參數為破裂起始點和破裂速度。由這兩個參數定義破裂傳播的規律。對于內陸地震，破裂規律主要和活斷層的特征相關［124］，其中俯沖型地震更多是參考歷史地震的破裂規律。

Recipe方法預測強地面運動的流程圖和震源參數設置步驟分別如圖1-2所示。由①②求得特定地震的震源模型；由②③來估計格林函數。利用震源模型和格林函數預測地震動，再由歷史地震數據驗證模擬地震動的合理性［125－130］。Recipe方法已經在多次歷史地震上成功驗證［131－132］，并且提出了多尺度空間異性的方法來驗證更大型的地震（≥M8.0）［133－134］。對于歷史記錄缺乏的區域，由于震源模型難以預估，采用不確定性來刻畫未知震源參數［135］，隨著方法不斷修正完善，Recipe 方法適用性會更強。不過Recipe 方法仍然存在缺陷如：凹凸體的數量無法確定，在沒有地震記錄的地區，震源模型參數如何設置，地層速度結構不夠精細等。

圖1 震源參數設置步驟Fig.1 The steps of setting source parameters

Recipe 方法利用了豐富的觀測資料和地層數據來刻畫震源參數，由此來模擬未來可能發生的地震動。Recipe 方法充分利用地震記錄豐富的優勢，通過大量地震數據反演斷層震源信息，在應用到未來可能發生的地震上。但是該方法在地震記錄較少的地方應用困難，而且該方法對斷層內部參數的設定仍需要完善。

圖2 Recipe方法預測地震動流程圖Fig.2 The flow chart of ground motion prediction by Recipe method

4 物理地震危險性分析中地震動數值模擬方法的展望

過去地震動的模擬或預測在面對地震觀測記錄稀少、對地層速度結構的認識有限和計算能力匱乏時，非常難以實現基于物理的地震動模擬和預測。想要實現大范圍的多線程的計算從計算上和理論上都是比較困難的。地震動數值模擬在近20 年來，發展極為迅速，不斷更新的計算模型和優越的計算條件都為應對復雜問題提供了基礎。未來，地震動數值模擬會是研究地震危險性的重要手段，甚至朝向取代經驗衰減關系的目標發展。雖然利用數值模擬在計算地震危險性上存在許多的不足與局限，但不否認其發展的潛力，在震害防御、抗震設防和工程結構抗震中都會有廣闊的發展前景。基于上述介紹，未來數值模擬在地震危險性分析中的研究及發展趨勢有以下幾個方面：

（1）將非彈性響應嵌入到對應計算程序中；構建更加精細的三維速度模型；更好地確定震源參數。在彈性介質中，如何將非彈性非線性的地震動波動現象高效和集中地模擬出來，是即將面對的一個難題。此外，還需要對三維速度結構（特別是淺表地層）做更深的探索，比如采用密集節點列陣、光纖傳感器或者是激光雷達的方式來完善地層結構模型，通過這些措施增加對斷層系統的理解，提高計算頻率等。對于基于凹凸體的地震動預測方法，如何更好的確定內在斷層參數（凹凸體位置和數量等參數），也是需要面臨的重大挑戰，這將會是未來研究的一個方向。

（2）與經驗衰減關系的博弈?；跀抵的M的地震危險性分析還處在一個研究、探索、驗證和完善的過程。那么不可避免地會與經驗衰減關系作比較，最近的研究表明數值模擬的方法是可行的。那么到底是用經驗關系還是地震動數值模擬的結果作為預測參考？本人認為目前還是以經驗關系為主，地震動數值模擬的方式作為修正。因為地震動數值模擬的方式尚處在發展和探索階段，并且研究區域較小，如若要擴大應用范圍、實際應用到工程之中需大量時間和基礎工作。因此，在未來，基于地震動數值模擬的地震危險性分析會成為計算地震危險性的常規手段，只不過還需不斷地發展、完善和論證。

（3）對我國的啟示。我國地處環太平洋地震帶與歐亞地震帶之間，地震活動十分活躍，同時也是全球地震災害較為嚴重的國家。隨著西部戰略和海洋強國戰略的實施，越來越多大型建筑都在山區及復雜地形環境中開展建設，那么抗震設計可能難以適應這些特殊地區的工程設防要求，而地震動數值模擬能夠很好地彌補這一缺陷，通過震源模型和速度結構模擬出地震波的傳播效應，從而為特殊區域提供地震危險性分析和提供設防依據。因此，地震動數值模擬會是研究區域地震危險性的重要手段之一。

數值模擬在物理地震危險性分析上的發展還遠沒有結束。作為涉及計算科學和地震科學等多學科、具有重要現實意義和廣泛應用的研究領域，研究地震危險性將會成為地震動數值模擬的重要研究方向。