茅遠(yuǎn)哲 呂國(guó)軍 孟立朋

1 河北省地震局，石家莊市槐中路262號(hào)，050021 2 河北省地震動(dòng)力學(xué)重點(diǎn)實(shí)驗(yàn)室，河北省三河市學(xué)院街465號(hào)，065201

地震的孕育和發(fā)生是孕震體應(yīng)力應(yīng)變能不斷積累、進(jìn)入臨界狀態(tài)并最終失穩(wěn)的力學(xué)過程，地殼巖石破裂過程中的應(yīng)力分布、應(yīng)力擾動(dòng)與斷裂構(gòu)造明顯相關(guān)[1]。因此，要探討強(qiáng)震的遷移規(guī)律并開展地震危險(xiǎn)性分析，最根本的途徑是研究斷裂帶的應(yīng)力狀態(tài)及其動(dòng)態(tài)演化過程[2]。

本文首先基于邯鄲地區(qū)活動(dòng)斷裂地表、淺層、中層、深層探測(cè)和小震重新定位等地球物理探測(cè)結(jié)果以及前人研究成果，建立該區(qū)域三維粘彈性有限元模擬的斷層模型；然后根據(jù)歷史地震和古地震研究結(jié)果，加入運(yùn)動(dòng)學(xué)特征邊界條件，模擬地震發(fā)生在目標(biāo)區(qū)活動(dòng)斷裂上引起的應(yīng)力應(yīng)變場(chǎng)變化；最后計(jì)算1830年磁縣地震引起的周圍斷裂面和滑動(dòng)方向上產(chǎn)生的庫(kù)侖破裂應(yīng)力變化，給出研究區(qū)活動(dòng)斷裂地震危險(xiǎn)性綜合評(píng)價(jià)結(jié)果。

1 計(jì)算方法與模型

1.1 目標(biāo)區(qū)內(nèi)主要斷裂分布

邯鄲地區(qū)主要斷裂呈北北東向和北西西向分布(圖1，藍(lán)色線段表示斷裂，紅色實(shí)心圓為磁縣地震)，北北東向的斷裂自西向東為紫山西斷裂、太行山山前斷裂和邯東斷裂(邯鄲隱伏斷裂)，北西西向斷裂自北向南為永年斷裂和磁縣斷裂。這些斷裂近似呈兩兩正交，控制著研究區(qū)內(nèi)的構(gòu)造活動(dòng)。研究區(qū)主要活斷層探測(cè)結(jié)果見表1。

圖1 目標(biāo)區(qū)主要構(gòu)造背景

表1 邯鄲地區(qū)主要斷裂要素

1.2 幾何模型的建立及介質(zhì)參數(shù)確定

1.2.1 幾何模型

根據(jù)地表和深部斷裂展布資料，建立目標(biāo)區(qū)幾何模型(圖2)。為防止邊界出現(xiàn)突變，要求模型區(qū)域比目標(biāo)區(qū)稍大，模型范圍以紫山西斷裂、太行山山前斷裂、邯東斷裂(邯鄲隱伏斷裂)、永年斷裂和磁縣斷裂為界，每層劃分為6個(gè)塊體。因目標(biāo)區(qū)內(nèi)地殼厚度起伏不大，設(shè)計(jì)模型每層水平展布，模型各層厚度見表2。

圖2 根據(jù)目標(biāo)區(qū)斷裂分布得到的平面幾何模型

表2 模型中不同深度塊體的橫波波速

由淺中層地震勘探獲得了主要斷裂的傾角，但在地殼深部這些主要斷裂的展布情況不明，因此在模型中將主要斷裂均設(shè)為近直立斷裂。由于目標(biāo)區(qū)主要在模型中心部分，因此在模型中部對(duì)網(wǎng)格進(jìn)行加密劃分，分成10 563個(gè)節(jié)點(diǎn)、8 736個(gè)單元(圖3，圖中藍(lán)色線段表示斷裂)。

1.2.2 模型介質(zhì)楊氏模量

利用有限元程序TEKTON進(jìn)行模擬計(jì)算[3]，該程序使用的計(jì)算方法為根據(jù)巖石上施加的應(yīng)力控制巖石破裂的機(jī)制，將莫爾圓實(shí)際半徑(O1-O3)/2與破裂時(shí)半徑之間的比值定義巖石極限破裂值Pf，可表示為：

式中，σ1與σ3分別為被施加的最大主應(yīng)力與最小主應(yīng)力。計(jì)算需要確定模型介質(zhì)的楊氏模量和粘滯性系數(shù)等參數(shù)。

利用深部探測(cè)結(jié)果來計(jì)算模型內(nèi)各不同深度塊體的楊氏模量。由于寬頻帶地震臺(tái)陣探測(cè)范圍比模型區(qū)域要小，主要在東西向跨太行山山前斷裂帶和邯東斷裂(邯鄲隱伏斷裂)、南北向跨永年斷裂和磁縣斷裂的區(qū)域內(nèi)布設(shè)；同時(shí)模型的主要關(guān)注區(qū)域在模型的中心部分。因此，用圖3中A1、A2、A5、A6塊體的波速代表整個(gè)塊體的波速(有一定程度的簡(jiǎn)化)。根據(jù)寬頻帶臺(tái)陣深部探測(cè)結(jié)果獲取不同深度模型塊體的橫波波速，見表2。

圖3 有限元模型結(jié)構(gòu)和單元?jiǎng)澐?/p>

地震波速與介質(zhì)彈性常數(shù)的關(guān)系為：

(1)

E=2μ(1+ν)

(2)

式中，E為楊氏模量，VS為S波速度，ρ為巖石密度，μ為剪切模量，ν為泊松比。根據(jù)地殼巖石的常用量，ν設(shè)定為0.25；根據(jù)文獻(xiàn)[4]，ρ取2.7×109g/m3。由式(1)可用VS和ρ求得μ，再根據(jù)式(2)求得E，由此得到有限元模型中不同深度塊體的楊氏模量(表3)。

表3 有限元模型不同深度塊體的楊氏模量

1.2.3 模型介質(zhì)有效粘滯性系數(shù)

周永勝等[5]在總結(jié)近30 a來高溫高壓流變實(shí)驗(yàn)資料的基礎(chǔ)上，應(yīng)用流變數(shù)據(jù)結(jié)合地震震源深度分布，對(duì)華北地殼流變性質(zhì)進(jìn)行了研究，得到華北地區(qū)巖石圈強(qiáng)度剖面模型中巖石圈各層的流變強(qiáng)度數(shù)據(jù)。本文根據(jù)其研究成果得出30 km深度范圍內(nèi)的有效粘滯性系數(shù)(表4)。根據(jù)GPS觀測(cè)的中國(guó)大陸地殼水平位移場(chǎng)計(jì)算出的最大剪應(yīng)變速率值，目標(biāo)區(qū)的剪切應(yīng)變率在1.268×10-15/s左右，但對(duì)于模型中30～80 km深度范圍，由于不確定其應(yīng)變速率是否與地表一致，因此應(yīng)變速率采用平均值1×10-15/s。

表4 有限元模型不同深度的有效粘滯性系數(shù)

1.2.4 模型邊界條件

影響目標(biāo)區(qū)動(dòng)力學(xué)特征的主要因素為巖石圈介質(zhì)特性和運(yùn)動(dòng)學(xué)特征，因此模擬時(shí)將目標(biāo)區(qū)的運(yùn)動(dòng)學(xué)特征作為邊界條件代入有限元模型進(jìn)行計(jì)算。

假設(shè)模型底部垂直于地表方向位移固定、地表自由，由于模型只到80 km深度，根據(jù)GPS位移場(chǎng)計(jì)算出模型4個(gè)邊界節(jié)點(diǎn)上的速度矢量(圖4)，其中目標(biāo)區(qū)內(nèi)部斷裂標(biāo)注較多，其GPS點(diǎn)位所示方向與區(qū)外一致，呈SE向。

圖4 華北地區(qū)GPS速度場(chǎng)

2 粘彈性有限元模型的建立和地震危險(xiǎn)性分析

利用得到的斷裂活動(dòng)模型、歷史強(qiáng)震數(shù)據(jù)以及區(qū)內(nèi)巖石圈介質(zhì)分布，可以建立粘彈性有限元模型，模擬地震發(fā)生在目標(biāo)區(qū)內(nèi)活動(dòng)斷裂上引起的應(yīng)力應(yīng)變場(chǎng)變化，計(jì)算活動(dòng)斷裂的斷裂面和滑動(dòng)方向上產(chǎn)生的庫(kù)侖破裂應(yīng)力變化，以分析目標(biāo)區(qū)未來強(qiáng)震危險(xiǎn)性。

1830年磁縣斷裂西段發(fā)生7.5級(jí)強(qiáng)震(表5)，利用斷層有限元計(jì)算中的劈節(jié)點(diǎn)技術(shù)[3]，在模擬的彈性步中按表5數(shù)據(jù)施加在磁縣斷裂上。

表5 1830年磁縣地震主要參數(shù)

研究區(qū)內(nèi)紫山西斷裂和磁縣斷裂中、西段現(xiàn)今時(shí)有小震發(fā)生，據(jù)此設(shè)計(jì)模型：磁縣強(qiáng)震發(fā)生，并且紫山西斷裂和磁縣斷裂中、西段小震不斷。

1830年磁縣地震距今已有190 a，而華北地區(qū)強(qiáng)震復(fù)發(fā)周期約在千年尺度，據(jù)此模擬的時(shí)間過程設(shè)計(jì)為：1)發(fā)震；2)震后幾天的調(diào)整；3)震后100 a的時(shí)間過程；4)震后1 000 a的時(shí)間過程，設(shè)計(jì)模擬的時(shí)間步為：1 d×2+2 a×2+10 a×10+100 a×3+1 000 a×3。

參照文獻(xiàn)[6-7]的方法，取μ′= 0.4。數(shù)值實(shí)驗(yàn)表明，改變此值對(duì)計(jì)算得到的庫(kù)侖破裂應(yīng)力變化的空間分布影響不大，但對(duì)其大小有一定的影響。將剪切應(yīng)力變化投影到假定的滑動(dòng)方向，與假定的滑動(dòng)方向一致取正，反之取負(fù)。假定的滑動(dòng)方向取自后續(xù)破裂事件的震源機(jī)制。

3 模擬結(jié)果與討論

根據(jù)以上模型設(shè)計(jì)，得到模型的模擬結(jié)果。首先根據(jù)時(shí)間尺度的不同，選取0 a、4 a、204 a、3 404 a的時(shí)間節(jié)點(diǎn)，提取時(shí)間步t0(彈性步)、t4(4 a)、t15(204 a)和t20(3 404 a)的位移場(chǎng)進(jìn)行分析(圖5、6)。可以看出，磁縣地震發(fā)生時(shí)(圖5)，由于存在水平和垂向位錯(cuò)，使得震中附近的水平位移場(chǎng)發(fā)生扭轉(zhuǎn)，相對(duì)于遠(yuǎn)離震中的位置，目標(biāo)區(qū)內(nèi)震中的位移量明顯大于其他地區(qū)。

圖5 磁縣地震時(shí)和震后4 a模型的位移場(chǎng)

震后4 a，目標(biāo)區(qū)內(nèi)不同位置位移量的差異減小，但區(qū)內(nèi)位移方向分布沒有明顯變化，說明強(qiáng)震引起的局部位移異常雖然在慢慢減弱，但是還未消除。震后204 a和3 404 a的位移圖像(圖6)表明，整個(gè)區(qū)域內(nèi)位移比較均勻，已經(jīng)看不到局部位移方向和大小的異常，說明磁縣地震的影響已經(jīng)被完全吸收。

圖6 磁縣地震震后204 a和3 404 a模型的位移場(chǎng)

圖7(單位MPa)為1830年磁縣地震引起的周圍斷裂面和滑動(dòng)方向上的庫(kù)侖破裂應(yīng)力變化。圖中藍(lán)色至綠色表示斷裂上庫(kù)侖應(yīng)力降低，地震危險(xiǎn)性降低；紅色至黃色表示斷裂上庫(kù)侖應(yīng)力增加，地震危險(xiǎn)性增加。由于紫山西斷裂和磁縣斷裂中、西段深部小震沿?cái)嗔逊植糩8]，確定為2條現(xiàn)今活動(dòng)斷裂，其庫(kù)侖應(yīng)力的變化對(duì)于地震危險(xiǎn)性評(píng)價(jià)至關(guān)重要。模型顯示，磁縣強(qiáng)震發(fā)生后，磁縣斷裂中、西段庫(kù)侖應(yīng)力加強(qiáng)，在震中附近位置庫(kù)侖應(yīng)力降低。對(duì)于紫山西斷裂，在斷裂南端，特別是與磁縣斷裂交接處庫(kù)侖應(yīng)力加強(qiáng)；在斷裂中段，庫(kù)侖應(yīng)力減弱；在斷裂中偏北的一小段斷裂上，庫(kù)侖應(yīng)力增強(qiáng)明顯；在斷裂北段，庫(kù)侖應(yīng)力增量值較小。說明在磁縣斷裂中、西段兩端與2條北東向斷裂交接的部位，地震危險(xiǎn)性值得關(guān)注。

圖7 磁縣地震在周圍活動(dòng)斷裂的斷裂面和滑動(dòng)方向上產(chǎn)生的庫(kù)侖破裂應(yīng)力變化

4 結(jié) 語(yǔ)

研究區(qū)內(nèi)，磁縣斷裂中、西段兩端與2條北東向斷裂(紫山西斷裂和太行山山前斷裂)交接處庫(kù)侖應(yīng)力有所加強(qiáng)，未來地震危險(xiǎn)性值得關(guān)注。

致謝：感謝中國(guó)地震局地質(zhì)研究所陶瑋副研究員為本研究提供巨大幫助。