茅遠哲 呂國軍 孟立朋
1 河北省地震局,石家莊市槐中路262號,050021 2 河北省地震動力學重點實驗室,河北省三河市學院街465號,065201
地震的孕育和發生是孕震體應力應變能不斷積累、進入臨界狀態并最終失穩的力學過程,地殼巖石破裂過程中的應力分布、應力擾動與斷裂構造明顯相關[1]。因此,要探討強震的遷移規律并開展地震危險性分析,最根本的途徑是研究斷裂帶的應力狀態及其動態演化過程[2]。
本文首先基于邯鄲地區活動斷裂地表、淺層、中層、深層探測和小震重新定位等地球物理探測結果以及前人研究成果,建立該區域三維粘彈性有限元模擬的斷層模型;然后根據歷史地震和古地震研究結果,加入運動學特征邊界條件,模擬地震發生在目標區活動斷裂上引起的應力應變場變化;最后計算1830年磁縣地震引起的周圍斷裂面和滑動方向上產生的庫侖破裂應力變化,給出研究區活動斷裂地震危險性綜合評價結果。
邯鄲地區主要斷裂呈北北東向和北西西向分布(圖1,藍色線段表示斷裂,紅色實心圓為磁縣地震),北北東向的斷裂自西向東為紫山西斷裂、太行山山前斷裂和邯東斷裂(邯鄲隱伏斷裂),北西西向斷裂自北向南為永年斷裂和磁縣斷裂。這些斷裂近似呈兩兩正交,控制著研究區內的構造活動。研究區主要活斷層探測結果見表1。

圖1 目標區主要構造背景

表1 邯鄲地區主要斷裂要素
1.2.1 幾何模型
根據地表和深部斷裂展布資料,建立目標區幾何模型(圖2)。為防止邊界出現突變,要求模型區域比目標區稍大,模型范圍以紫山西斷裂、太行山山前斷裂、邯東斷裂(邯鄲隱伏斷裂)、永年斷裂和磁縣斷裂為界,每層劃分為6個塊體。因目標區內地殼厚度起伏不大,設計模型每層水平展布,模型各層厚度見表2。

圖2 根據目標區斷裂分布得到的平面幾何模型

表2 模型中不同深度塊體的橫波波速
由淺中層地震勘探獲得了主要斷裂的傾角,但在地殼深部這些主要斷裂的展布情況不明,因此在模型中將主要斷裂均設為近直立斷裂。由于目標區主要在模型中心部分,因此在模型中部對網格進行加密劃分,分成10 563個節點、8 736個單元(圖3,圖中藍色線段表示斷裂)。
1.2.2 模型介質楊氏模量
利用有限元程序TEKTON進行模擬計算[3],該程序使用的計算方法為根據巖石上施加的應力控制巖石破裂的機制,將莫爾圓實際半徑(O1-O3)/2與破裂時半徑之間的比值定義巖石極限破裂值Pf,可表示為:
式中,σ1與σ3分別為被施加的最大主應力與最小主應力。計算需要確定模型介質的楊氏模量和粘滯性系數等參數。
利用深部探測結果來計算模型內各不同深度塊體的楊氏模量。由于寬頻帶地震臺陣探測范圍比模型區域要小,主要在東西向跨太行山山前斷裂帶和邯東斷裂(邯鄲隱伏斷裂)、南北向跨永年斷裂和磁縣斷裂的區域內布設;同時模型的主要關注區域在模型的中心部分。因此,用圖3中A1、A2、A5、A6塊體的波速代表整個塊體的波速(有一定程度的簡化)。根據寬頻帶臺陣深部探測結果獲取不同深度模型塊體的橫波波速,見表2。

圖3 有限元模型結構和單元劃分
地震波速與介質彈性常數的關系為:
(1)
E=2μ(1+ν)
(2)
式中,E為楊氏模量,VS為S波速度,ρ為巖石密度,μ為剪切模量,ν為泊松比。根據地殼巖石的常用量,ν設定為0.25;根據文獻[4],ρ取2.7×109g/m3。由式(1)可用VS和ρ求得μ,再根據式(2)求得E,由此得到有限元模型中不同深度塊體的楊氏模量(表3)。

表3 有限元模型不同深度塊體的楊氏模量
1.2.3 模型介質有效粘滯性系數
周永勝等[5]在總結近30 a來高溫高壓流變實驗資料的基礎上,應用流變數據結合地震震源深度分布,對華北地殼流變性質進行了研究,得到華北地區巖石圈強度剖面模型中巖石圈各層的流變強度數據。本文根據其研究成果得出30 km深度范圍內的有效粘滯性系數(表4)。根據GPS觀測的中國大陸地殼水平位移場計算出的最大剪應變速率值,目標區的剪切應變率在1.268×10-15/s左右,但對于模型中30~80 km深度范圍,由于不確定其應變速率是否與地表一致,因此應變速率采用平均值1×10-15/s。

表4 有限元模型不同深度的有效粘滯性系數
1.2.4 模型邊界條件
影響目標區動力學特征的主要因素為巖石圈介質特性和運動學特征,因此模擬時將目標區的運動學特征作為邊界條件代入有限元模型進行計算。
假設模型底部垂直于地表方向位移固定、地表自由,由于模型只到80 km深度,根據GPS位移場計算出模型4個邊界節點上的速度矢量(圖4),其中目標區內部斷裂標注較多,其GPS點位所示方向與區外一致,呈SE向。

圖4 華北地區GPS速度場
利用得到的斷裂活動模型、歷史強震數據以及區內巖石圈介質分布,可以建立粘彈性有限元模型,模擬地震發生在目標區內活動斷裂上引起的應力應變場變化,計算活動斷裂的斷裂面和滑動方向上產生的庫侖破裂應力變化,以分析目標區未來強震危險性。
1830年磁縣斷裂西段發生7.5級強震(表5),利用斷層有限元計算中的劈節點技術[3],在模擬的彈性步中按表5數據施加在磁縣斷裂上。

表5 1830年磁縣地震主要參數
研究區內紫山西斷裂和磁縣斷裂中、西段現今時有小震發生,據此設計模型:磁縣強震發生,并且紫山西斷裂和磁縣斷裂中、西段小震不斷。
1830年磁縣地震距今已有190 a,而華北地區強震復發周期約在千年尺度,據此模擬的時間過程設計為:1)發震;2)震后幾天的調整;3)震后100 a的時間過程;4)震后1 000 a的時間過程,設計模擬的時間步為:1 d×2+2 a×2+10 a×10+100 a×3+1 000 a×3。
參照文獻[6-7]的方法,取μ′= 0.4。數值實驗表明,改變此值對計算得到的庫侖破裂應力變化的空間分布影響不大,但對其大小有一定的影響。將剪切應力變化投影到假定的滑動方向,與假定的滑動方向一致取正,反之取負。假定的滑動方向取自后續破裂事件的震源機制。
根據以上模型設計,得到模型的模擬結果。首先根據時間尺度的不同,選取0 a、4 a、204 a、3 404 a的時間節點,提取時間步t0(彈性步)、t4(4 a)、t15(204 a)和t20(3 404 a)的位移場進行分析(圖5、6)。可以看出,磁縣地震發生時(圖5),由于存在水平和垂向位錯,使得震中附近的水平位移場發生扭轉,相對于遠離震中的位置,目標區內震中的位移量明顯大于其他地區。

圖5 磁縣地震時和震后4 a模型的位移場
震后4 a,目標區內不同位置位移量的差異減小,但區內位移方向分布沒有明顯變化,說明強震引起的局部位移異常雖然在慢慢減弱,但是還未消除。震后204 a和3 404 a的位移圖像(圖6)表明,整個區域內位移比較均勻,已經看不到局部位移方向和大小的異常,說明磁縣地震的影響已經被完全吸收。

圖6 磁縣地震震后204 a和3 404 a模型的位移場
圖7(單位MPa)為1830年磁縣地震引起的周圍斷裂面和滑動方向上的庫侖破裂應力變化。圖中藍色至綠色表示斷裂上庫侖應力降低,地震危險性降低;紅色至黃色表示斷裂上庫侖應力增加,地震危險性增加。由于紫山西斷裂和磁縣斷裂中、西段深部小震沿斷裂分布[8],確定為2條現今活動斷裂,其庫侖應力的變化對于地震危險性評價至關重要。模型顯示,磁縣強震發生后,磁縣斷裂中、西段庫侖應力加強,在震中附近位置庫侖應力降低。對于紫山西斷裂,在斷裂南端,特別是與磁縣斷裂交接處庫侖應力加強;在斷裂中段,庫侖應力減弱;在斷裂中偏北的一小段斷裂上,庫侖應力增強明顯;在斷裂北段,庫侖應力增量值較小。說明在磁縣斷裂中、西段兩端與2條北東向斷裂交接的部位,地震危險性值得關注。

圖7 磁縣地震在周圍活動斷裂的斷裂面和滑動方向上產生的庫侖破裂應力變化
研究區內,磁縣斷裂中、西段兩端與2條北東向斷裂(紫山西斷裂和太行山山前斷裂)交接處庫侖應力有所加強,未來地震危險性值得關注。
致謝:感謝中國地震局地質研究所陶瑋副研究員為本研究提供巨大幫助。