三峽水庫(kù)蓄水后仙女山斷裂北段數(shù)值模擬及其孕震機(jī)理研究

吳建超　陳蜀俊　陳俊華　蔡永建

摘要：基于線彈性理論有限元方法，建立三峽庫(kù)首區(qū)的三維有限元模型，模擬計(jì)算了庫(kù)首區(qū)尤其是仙女山斷裂北段在蓄水后的全位移場(chǎng)和構(gòu)造應(yīng)力場(chǎng)變化。結(jié)果表明：三峽水庫(kù)175m水位蓄水，造成仙女山斷裂和九畹溪斷裂北段處于全位移場(chǎng)歸州低值區(qū)和新灘高值區(qū)的梯度帶上：同時(shí)，仙女山斷裂和九畹溪斷裂北段深部3組主應(yīng)力值較蓄水前陡增。這樣的地段全位移差值變化明顯，3組主應(yīng)力值增高，剪應(yīng)力和剪應(yīng)變能易于集中，容易加速仙女山斷裂和九畹溪斷裂北端失穩(wěn)而誘發(fā)水庫(kù)地震。

關(guān)鍵詞：三峽水庫(kù)；仙女山斷裂；數(shù)值模擬；孕震機(jī)理；水庫(kù)地震

中圖分類(lèi)號(hào)：P315.7 文獻(xiàn)標(biāo)識(shí)碼：A 文章編號(hào)：1000-0666（2016）02-0218-06

0 引言

三峽水庫(kù)自2003年6月蓄水以來(lái)，在庫(kù)區(qū)周邊發(fā)生了數(shù)千次地震，絕大多數(shù)為M≤2.0非構(gòu)造成因的誘發(fā)地震（車(chē)用太等，2009）。然而，伴隨著三峽水庫(kù)175m試驗(yàn)性蓄水進(jìn)程，2008年11月22日、2014年3月27日和3月30日，仙女山斷裂北段附近分別發(fā)生了M_S4.1、M_S4.2和M_S4.5地震，引發(fā)了社會(huì)各界對(duì)于三峽水庫(kù)是否誘發(fā)了地震的討論。

筆者參加了2014年3月30日秭歸M_S4.5地震的應(yīng)急調(diào)查工作，該次地震震中烈度為Ⅴ度，有感范圍較大。在郭家壩鎮(zhèn)頭道河村、郭家壩村等地出現(xiàn)Ⅵ度烈度異常區(qū)，且造成了一定程度的房屋破壞和財(cái)產(chǎn)損失。根據(jù)地震監(jiān)測(cè)記錄波形分析，認(rèn)為此次地震為構(gòu)造型水庫(kù)地震，與仙女山微地塊的活動(dòng)有關(guān)。

值得注意的是，蓄水后三峽庫(kù)區(qū)的地震活動(dòng)多發(fā)生在仙女山斷裂和九畹溪斷裂的圍限區(qū)附近，屬于這兩條斷裂的北端部。斷層端部構(gòu)造應(yīng)力場(chǎng)及孕震機(jī)理一直是國(guó)內(nèi)外研究的熱點(diǎn)。斷層的末端擴(kuò)展受阻，往往會(huì)在斷層的端部產(chǎn)生應(yīng)力集中而引發(fā)地震（Obara et al，2009；孔禮健等，2009）。而三峽水庫(kù)蓄水后的形變觀測(cè)結(jié)果顯示，仙女山斷裂北段地塊的傾斜場(chǎng)發(fā)生了改變，由蓄水前的正斷層運(yùn)動(dòng)轉(zhuǎn)為蓄水后的逆斷運(yùn)動(dòng)性質(zhì)（張燕等，2005；邢燦飛等，2006）。因此，綜合考慮蓄水效應(yīng)，研究蓄水后仙女山斷裂北段構(gòu)造應(yīng)力場(chǎng)及其孕震機(jī)理具有重要意義。

本文應(yīng)用三維有限元方法模擬計(jì)算了三峽水庫(kù)蓄水對(duì)仙女山斷裂北段構(gòu)造應(yīng)力場(chǎng)的影響，并對(duì)其孕震機(jī)理進(jìn)行了分析和討論，以期對(duì)三峽水庫(kù)誘發(fā)地震機(jī)理作進(jìn)一步探索。

1 地震活動(dòng)概況

三峽水庫(kù)自2003年6月蓄水以來(lái)，仙女山斷裂北段附近共發(fā)生了41次M_S≥2.0地震，其中M_S4.0～4.9地震3次，M_S3.0～3.9地震3次，M_S2.0～2.9地震35次，且大部分位于長(zhǎng)江南岸，呈密集面狀分布特征，如圖1所示，地震目錄取自湖北省地震局臺(tái)網(wǎng)中心。

車(chē)用太等（2009）認(rèn)為，2008年11月22日胡家坪M_S4.1地震是在水庫(kù)水體荷載與庫(kù)水下滲的共同作用下沿仙女山斷裂發(fā)生的構(gòu)造型水庫(kù)誘發(fā)地震；吳建超等（2012）認(rèn)為胡家坪M_S4.1地震發(fā)生在等效應(yīng)力高值區(qū)與低值區(qū)的梯度帶上，為等效應(yīng)力變化最為劇烈、最易產(chǎn)生應(yīng)變能積累和應(yīng)力集中的區(qū)域；陳俊華等（2013）從震群活動(dòng)時(shí)序特征、震源體大小、震源機(jī)制及地震波形的頻譜特征等方面進(jìn)行研究，認(rèn)為2012年10月31日屈原鎮(zhèn)M_S3.2地震為水庫(kù)觸發(fā)型構(gòu)造地震；徐長(zhǎng)朋（2012）對(duì)三峽庫(kù)區(qū)119次小震的震源機(jī)制解進(jìn)行研究，認(rèn)為仙女山斷裂北段地震震源機(jī)制以走滑為主。

上述研究多針對(duì)蓄水所誘發(fā)的地震，未對(duì)蓄水后仙女山斷裂上的位移場(chǎng)和深部應(yīng)力場(chǎng)變化進(jìn)行綜合研究。因此，有必要通過(guò)對(duì)蓄水后仙女山斷裂北段的孕震動(dòng)力學(xué)環(huán)境進(jìn)行數(shù)值模擬，從上地殼深部探究其孕震機(jī)理。

2 研究區(qū)三維有限元模型

本文構(gòu)建的三維有限元模型范圍為110km×85km，深度取至上地殼底部埋深12km處。魏貴春等（2013）利用雙差定位法對(duì)三峽數(shù)字遙測(cè)地震臺(tái)網(wǎng)記錄到的秭歸附近地區(qū)73個(gè)M_S≥1.0的地震進(jìn)行了重新定位，結(jié)果顯示，絕大多數(shù)地震事件的震源深度均在7km左右。這一深度剛好位于本文三維有限元模型的中部，計(jì)算過(guò)程中避免了受邊界條件的影響，結(jié)果相對(duì)準(zhǔn)確。因此，本文構(gòu)建的三維有限元模型是合理與可靠的。圖2為考慮地形地貌、斷裂和構(gòu)造層建立的三維模型，模型以（110°00′00″N，31°56′35″E）為坐標(biāo)原點(diǎn)，X軸正方向?yàn)镾E33°，Y軸正方向?yàn)镹E33°，z軸垂直于水平面，向上為正。

模型按上地殼巖石層的物性參數(shù)可劃分為沉積蓋層和基底構(gòu)造層。其中，沉積蓋層厚約7km，主要為花崗巖、灰?guī)r和砂巖；基底構(gòu)造層層厚約5km，主要為花崗巖。依據(jù)陳蜀?。?005）研究成果，斷裂帶的參數(shù)分別取為：彈性模量E=1.07×10¹⁰Pa，泊松比v=0.24，密度ρ=2570kg·m^-3。有限元網(wǎng)格劃分以滿足計(jì)算精度為要求，對(duì)仙女山斷裂北段附近的網(wǎng)格再做細(xì)化。本次計(jì)算模型單元數(shù)共計(jì)18862個(gè)、節(jié)點(diǎn)數(shù)共計(jì)22498個(gè)。

參照陳蜀俊等（2005，2006）的研究成果和三峽地區(qū)地震加密觀測(cè)資料（馬文濤等，2010），本文選用的加載方案為：主壓應(yīng)力方向?yàn)镹E-SW向，大小為790MPa，主張應(yīng)力方向?yàn)镹W-SE向，大小為540MPa。有限元模型采用三維線彈性模型，除作了約束的部位外，每個(gè)節(jié)點(diǎn)有3個(gè)自由度。

3 計(jì)算結(jié)果分析與討論

三峽水庫(kù)蓄水提升了長(zhǎng)江上游水位，并淹沒(méi)了長(zhǎng)江部分支流河段。對(duì)庫(kù)水影響范圍內(nèi)的水文地質(zhì)條件與巖石物性參數(shù)進(jìn)行研究，是蓄水背景下構(gòu)造應(yīng)力場(chǎng)數(shù)值模擬的基礎(chǔ)。根據(jù)車(chē)用太等（1998）、劉遠(yuǎn)征等（2010）對(duì)水庫(kù)誘發(fā)地震成因機(jī)理的研究成果，本文著重考慮了三峽水庫(kù)蓄水后，水體荷載作用、庫(kù)水下滲作用和孔隙水壓作用的影響，并計(jì)算得到了蓄水后三峽庫(kù)首區(qū)的全位移場(chǎng)和構(gòu)造應(yīng)力場(chǎng)。

3.1 全位移場(chǎng)

為更具體地分析三峽水庫(kù)蓄水對(duì)仙女山斷裂北段全位移場(chǎng)的影響，將175m水位與蓄水前全位移值相減，得到了三峽水庫(kù)175m水位相對(duì)于蓄水前的全位移差值圖（圖3），與杜瑞林等（2004）給出的蓄水后GPS形變場(chǎng)觀測(cè)結(jié)果較為一致。

蓄水造成的垂直形變遠(yuǎn)大于水平形變。因此，垂直位移對(duì)全位移作出的貢獻(xiàn)最大。由圖3可知，175m水位相對(duì)于蓄水前的全位移差值分布不均勻，出現(xiàn)了歸州和官渡口等幾個(gè)低值區(qū)，最大值出現(xiàn)在巴東附近，而最小值出現(xiàn)在歸州一帶。極值區(qū)大部分位于長(zhǎng)江南岸，緊鄰極值點(diǎn)的等值線梯度大，其后逐漸變緩，極值區(qū)等值線長(zhǎng)軸方向與斷裂帶走向基本一致。值得注意的是，仙女山斷裂和九畹溪斷裂北段剛好處于歸州低值區(qū)和新灘高值區(qū)的梯度帶上。這樣的地段全位移變化明顯，易于剪應(yīng)力和垂直剪應(yīng)變能的集中，斷裂容易失穩(wěn)錯(cuò)動(dòng)而引發(fā)地震。

將模擬結(jié)果與觀測(cè)結(jié)果進(jìn)行對(duì)比分析可知，數(shù)值模擬獲得全位移場(chǎng)與GPS觀測(cè)形變場(chǎng)所反映的研究區(qū)垂直形變特征及趨勢(shì)基本吻合，全位移場(chǎng)的極低值區(qū)在歸州-新灘-帶出現(xiàn)（杜瑞林等，2004）。雖然上述模擬結(jié)果與GPS觀測(cè)結(jié)果所反映的運(yùn)動(dòng)趨勢(shì)整體上基本一致，但在局部區(qū)域還是存在差異，例如模型的邊界等。這可能是由于有限元計(jì)算中對(duì)于模型邊界的加載和圍限作用所致，而其對(duì)處于模型中部的仙女山斷裂附近則影響較小。因此，本文關(guān)心的仙女山斷裂北段的數(shù)值模擬結(jié)果是可信的。

3.2 構(gòu)造應(yīng)力場(chǎng)

蓄水前后構(gòu)造應(yīng)力場(chǎng)的變化，直接反映了由蓄水引起的構(gòu)造應(yīng)力的變化，是蓄水引起庫(kù)首區(qū)孕震環(huán)境變化的有力證據(jù)之一。為了詳細(xì)了解三峽水庫(kù)在蓄水進(jìn)程中產(chǎn)生的主應(yīng)力場(chǎng)變化，尤其是對(duì)仙女山斷裂深部應(yīng)力變化進(jìn)行詳細(xì)的剖析和研究，本文對(duì)模擬結(jié)果做了剖面輸出。以坐標(biāo)點(diǎn)（39400，0，-8000）作為切面點(diǎn)，以向量{-1，1，0}作為切面的法向向量，切出了跨仙女山斷裂和九畹溪斷裂的剖面線（圖1），蓄水前后跨仙女山斷裂和九畹溪斷裂的最大主應(yīng)力結(jié)果如圖4所示。

由圖4可以看出，伴隨著蓄水進(jìn)程，跨仙女山斷裂北段最大主應(yīng)力較蓄水前出現(xiàn)明顯變化。蓄水前，深部最大主應(yīng)力在郭家壩至新灘一帶沒(méi)有明顯變化，呈舒緩波狀分布，應(yīng)力值自上而下緩慢遞增。而蓄水至135m水位后，最大主應(yīng)力等值線在郭家壩和新灘間，即仙女山斷裂和九畹溪斷裂深部出現(xiàn)了較大變化，具體表現(xiàn)為兩組斷裂帶深部的最大主應(yīng)力陡增，呈“V”字形，斷裂帶內(nèi)部最大主應(yīng)力均大于斷裂帶外部的應(yīng)力大小。在175m水位下模擬計(jì)算得到的最大主應(yīng)力云圖基本同135m水位下的類(lèi)似，只是在模型頂部有應(yīng)力值的變化。

本文的研究重點(diǎn)是仙女山斷裂北段深部構(gòu)造應(yīng)力場(chǎng)的變化，考慮到該斷裂北段地震的震源深度多為7km左右。因此，取Z=-7km作為輸出層，得到了蓄水至175m水位后，庫(kù)首區(qū)7km深度處跨仙女山斷裂的主應(yīng)力場(chǎng)的結(jié)果。從圖5可以看出，175m水位蓄水后，在仙女山斷裂和九畹溪斷裂帶上均出現(xiàn)了3組主應(yīng)力大小陡增的現(xiàn)象，且均大于斷裂帶兩側(cè)的主應(yīng)力值。而兩條斷裂帶中間所挾持的區(qū)域則明顯地出現(xiàn)了主應(yīng)力減小，總體呈現(xiàn)“M”型。

上述數(shù)值模擬結(jié)果集中表明，三峽水庫(kù)蓄水給包括仙女山斷裂和九畹溪斷裂的整個(gè)庫(kù)首區(qū)的應(yīng)力狀態(tài)和位移場(chǎng)帶來(lái)了明顯變化。

4 應(yīng)力場(chǎng)變化對(duì)孕震環(huán)境的影響

大地構(gòu)造上，仙女山斷裂北段位于黃陵背斜的西南緣。黃陵背斜為一古老的剛性地質(zhì)體，在地質(zhì)歷史中長(zhǎng)期表現(xiàn)為隆起，對(duì)局部應(yīng)力作用產(chǎn)生強(qiáng)大的反作用力。因此，黃陵背斜周緣的隆起區(qū)存在垂直向上的作用力，從而導(dǎo)致上部地殼引張（高士鈞，1992）。

在三峽水庫(kù)蓄水進(jìn)程中，水體荷載會(huì)在庫(kù)底巖體內(nèi)產(chǎn)生彈性附加應(yīng)力的響應(yīng)，進(jìn)而將影響到周邊斷層穩(wěn)定性。在不同的構(gòu)造應(yīng)力環(huán)境中，附加應(yīng)力對(duì)斷層穩(wěn)定性所起的作用是截然不同的。在拉張應(yīng)力環(huán)境中，最大主應(yīng)力方向近于垂直，最小主應(yīng)力方向近于水平，水體荷載與初始構(gòu)造應(yīng)力場(chǎng)疊加，將會(huì)使莫爾圓半徑變大，增加逼近外包絡(luò)線的機(jī)會(huì)，進(jìn)而加劇斷層破裂而發(fā)生地震（周斌等，2010）。

該理論用于復(fù)雜的三維斷面并不簡(jiǎn)單，而用相近的辦法和思路，也可獲得相似的認(rèn)識(shí)。即蓄水導(dǎo)致斷裂帶上構(gòu)造應(yīng)力場(chǎng)發(fā)生變化，并疊加在原有的構(gòu)造應(yīng)力場(chǎng)上，形成了對(duì)構(gòu)造應(yīng)力場(chǎng)及孕震環(huán)境的擾動(dòng)和影響，將產(chǎn)生軟化、促滑、誘震作用，使原本相對(duì)穩(wěn)定的斷裂發(fā)生活動(dòng)或原本就存在的地震危險(xiǎn)趨勢(shì)增加，甚至使未來(lái)地震震源機(jī)制和時(shí)、空、強(qiáng)度特征發(fā)生改變。

仙女山斷裂北段主應(yīng)力的模擬計(jì)算結(jié)果表明，蓄水后斷裂帶上的應(yīng)力加載效應(yīng)顯著。一般而言，加載效應(yīng)有利于局部應(yīng)力積累，可能加速地震發(fā)生或提高局部地震危險(xiǎn)度。卸載效應(yīng)有利于區(qū)域應(yīng)力釋放，可能延遲地震發(fā)生或降低區(qū)域地震能量。

張秋文等（2005）研究認(rèn)為，幾個(gè)兆帕的應(yīng)力場(chǎng)變化，雖不足以直接導(dǎo)致地震，卻可以使地震發(fā)生的時(shí)間提前、發(fā)震概率提高。王紹晉等（2005）研究發(fā)現(xiàn)，瀾滄江漫灣水庫(kù)蓄水前3.0≤M_S≤4.0的歷史地震震中均在庫(kù)中心10km以外，震源機(jī)制解平均應(yīng)力場(chǎng)主壓應(yīng)力方位為SSE或近SN向，蓄水后震中向庫(kù)心匯集，震源機(jī)制解表明斷層面錯(cuò)動(dòng)具有較大的正傾滑動(dòng)分量，主應(yīng)力方位則變?yōu)榻麰W向，方位向東偏轉(zhuǎn)了50°以上。

綜上所述，三峽水庫(kù)175m水位蓄水后，仙女山斷裂北段淺部構(gòu)造應(yīng)力場(chǎng)的動(dòng)態(tài)變化、調(diào)整蘊(yùn)含著足以引起重視的地震動(dòng)力學(xué)意義。仙女山斷裂和九畹溪斷裂的北段之間全位移變化明顯，3組主應(yīng)力大小陡增，應(yīng)力加載效應(yīng)明顯，容易加速仙女山斷裂北端部失穩(wěn)而引發(fā)地震。

5 結(jié)論

本文模擬計(jì)算結(jié)果雖然在局部區(qū)域和模型邊界上與觀測(cè)結(jié)果存在一定的差異，但總體變化趨勢(shì)基本一致。仙女山斷裂北段深部構(gòu)造應(yīng)力場(chǎng)的變化，能夠綜合反映蓄水效應(yīng)，而構(gòu)造應(yīng)力場(chǎng)局部應(yīng)力狀態(tài)的調(diào)整又與地震的孕育存在直接聯(lián)系。因此，本文通過(guò)三維數(shù)值模擬研究蓄水后庫(kù)首區(qū)構(gòu)造應(yīng)力場(chǎng)的變化，將有助于探討蓄水后斷裂帶深部的孕震機(jī)理及其地震動(dòng)力學(xué)意義。

通過(guò)對(duì)計(jì)算結(jié)果的分析討論，我們獲得了以下主要認(rèn)識(shí)：

（1）利用三維有限元方法，計(jì)算得到了三峽庫(kù)區(qū)175m水位相對(duì)于蓄水前的全位移差值結(jié)果，與蓄水后GPS形變場(chǎng)觀測(cè)結(jié)果較為一致。沿歸州、新灘和周坪三地出現(xiàn)了全位移場(chǎng)低值區(qū)，反映了175m水位蓄水造成的上地殼表層全位移場(chǎng)變化的基本狀況。

（2）伴隨三峽水庫(kù)135m和175m水位的蓄水進(jìn)程，仙女山斷裂和九畹溪斷裂淺部最大主應(yīng)力較蓄水前出現(xiàn)明顯變化，表現(xiàn)為兩組斷裂帶深部的最大主應(yīng)力值陡增，且均大于斷裂帶兩側(cè)的數(shù)值。

（3）三峽水庫(kù)蓄水導(dǎo)致了仙女山斷裂北段淺部構(gòu)造應(yīng)力場(chǎng)和全位移場(chǎng)的動(dòng)態(tài)變化、調(diào)整。蓄水后，地震多發(fā)生于仙女山斷裂和九畹溪斷裂的北端附近，應(yīng)力易于在這兩組斷裂端部積累。

（4）仙女山斷裂北段在三峽地區(qū)的位置重要，175m水位蓄水以來(lái)，庫(kù)首區(qū)M≥2.0地震多發(fā)生于這個(gè)區(qū)域。綜合地層巖性、斷裂構(gòu)造、滲透條件等，運(yùn)用水-巖耦合數(shù)值模擬手段，探究該區(qū)域地震的成因機(jī)理是我們下一步研究的目標(biāo)。