三峽水庫蓄水后仙女山斷裂北段數值模擬及其孕震機理研究

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

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

關鍵詞：三峽水庫；仙女山斷裂；數值模擬；孕震機理；水庫地震

中圖分類號：P315.7 文獻標識碼：A 文章編號：1000-0666（2016）02-0218-06

0 引言

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

筆者參加了2014年3月30日秭歸M_S4.5地震的應急調查工作，該次地震震中烈度為Ⅴ度，有感范圍較大。在郭家壩鎮頭道河村、郭家壩村等地出現Ⅵ度烈度異常區，且造成了一定程度的房屋破壞和財產損失。根據地震監測記錄波形分析，認為此次地震為構造型水庫地震，與仙女山微地塊的活動有關。

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

本文應用三維有限元方法模擬計算了三峽水庫蓄水對仙女山斷裂北段構造應力場的影響，并對其孕震機理進行了分析和討論，以期對三峽水庫誘發地震機理作進一步探索。

1 地震活動概況

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

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

上述研究多針對蓄水所誘發的地震，未對蓄水后仙女山斷裂上的位移場和深部應力場變化進行綜合研究。因此，有必要通過對蓄水后仙女山斷裂北段的孕震動力學環境進行數值模擬，從上地殼深部探究其孕震機理。

2 研究區三維有限元模型

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

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

參照陳蜀俊等（2005，2006）的研究成果和三峽地區地震加密觀測資料（馬文濤等，2010），本文選用的加載方案為：主壓應力方向為NE-SW向，大小為790MPa，主張應力方向為NW-SE向，大小為540MPa。有限元模型采用三維線彈性模型，除作了約束的部位外，每個節點有3個自由度。

3 計算結果分析與討論

三峽水庫蓄水提升了長江上游水位，并淹沒了長江部分支流河段。對庫水影響范圍內的水文地質條件與巖石物性參數進行研究，是蓄水背景下構造應力場數值模擬的基礎。根據車用太等（1998）、劉遠征等（2010）對水庫誘發地震成因機理的研究成果，本文著重考慮了三峽水庫蓄水后，水體荷載作用、庫水下滲作用和孔隙水壓作用的影響，并計算得到了蓄水后三峽庫首區的全位移場和構造應力場。

3.1 全位移場

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

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

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

3.2 構造應力場

蓄水前后構造應力場的變化，直接反映了由蓄水引起的構造應力的變化，是蓄水引起庫首區孕震環境變化的有力證據之一。為了詳細了解三峽水庫在蓄水進程中產生的主應力場變化，尤其是對仙女山斷裂深部應力變化進行詳細的剖析和研究，本文對模擬結果做了剖面輸出。以坐標點（39400，0，-8000）作為切面點，以向量{-1，1，0}作為切面的法向向量，切出了跨仙女山斷裂和九畹溪斷裂的剖面線（圖1），蓄水前后跨仙女山斷裂和九畹溪斷裂的最大主應力結果如圖4所示。

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

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

上述數值模擬結果集中表明，三峽水庫蓄水給包括仙女山斷裂和九畹溪斷裂的整個庫首區的應力狀態和位移場帶來了明顯變化。

4 應力場變化對孕震環境的影響

大地構造上，仙女山斷裂北段位于黃陵背斜的西南緣。黃陵背斜為一古老的剛性地質體，在地質歷史中長期表現為隆起，對局部應力作用產生強大的反作用力。因此，黃陵背斜周緣的隆起區存在垂直向上的作用力，從而導致上部地殼引張（高士鈞，1992）。

在三峽水庫蓄水進程中，水體荷載會在庫底巖體內產生彈性附加應力的響應，進而將影響到周邊斷層穩定性。在不同的構造應力環境中，附加應力對斷層穩定性所起的作用是截然不同的。在拉張應力環境中，最大主應力方向近于垂直，最小主應力方向近于水平，水體荷載與初始構造應力場疊加，將會使莫爾圓半徑變大，增加逼近外包絡線的機會，進而加劇斷層破裂而發生地震（周斌等，2010）。

該理論用于復雜的三維斷面并不簡單，而用相近的辦法和思路，也可獲得相似的認識。即蓄水導致斷裂帶上構造應力場發生變化，并疊加在原有的構造應力場上，形成了對構造應力場及孕震環境的擾動和影響，將產生軟化、促滑、誘震作用，使原本相對穩定的斷裂發生活動或原本就存在的地震危險趨勢增加，甚至使未來地震震源機制和時、空、強度特征發生改變。

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

張秋文等（2005）研究認為，幾個兆帕的應力場變化，雖不足以直接導致地震，卻可以使地震發生的時間提前、發震概率提高。王紹晉等（2005）研究發現，瀾滄江漫灣水庫蓄水前3.0≤M_S≤4.0的歷史地震震中均在庫中心10km以外，震源機制解平均應力場主壓應力方位為SSE或近SN向，蓄水后震中向庫心匯集，震源機制解表明斷層面錯動具有較大的正傾滑動分量，主應力方位則變為近EW向，方位向東偏轉了50°以上。

綜上所述，三峽水庫175m水位蓄水后，仙女山斷裂北段淺部構造應力場的動態變化、調整蘊含著足以引起重視的地震動力學意義。仙女山斷裂和九畹溪斷裂的北段之間全位移變化明顯，3組主應力大小陡增，應力加載效應明顯，容易加速仙女山斷裂北端部失穩而引發地震。

5 結論

本文模擬計算結果雖然在局部區域和模型邊界上與觀測結果存在一定的差異，但總體變化趨勢基本一致。仙女山斷裂北段深部構造應力場的變化，能夠綜合反映蓄水效應，而構造應力場局部應力狀態的調整又與地震的孕育存在直接聯系。因此，本文通過三維數值模擬研究蓄水后庫首區構造應力場的變化，將有助于探討蓄水后斷裂帶深部的孕震機理及其地震動力學意義。

通過對計算結果的分析討論，我們獲得了以下主要認識：

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

（2）伴隨三峽水庫135m和175m水位的蓄水進程，仙女山斷裂和九畹溪斷裂淺部最大主應力較蓄水前出現明顯變化，表現為兩組斷裂帶深部的最大主應力值陡增，且均大于斷裂帶兩側的數值。

（3）三峽水庫蓄水導致了仙女山斷裂北段淺部構造應力場和全位移場的動態變化、調整。蓄水后，地震多發生于仙女山斷裂和九畹溪斷裂的北端附近，應力易于在這兩組斷裂端部積累。

（4）仙女山斷裂北段在三峽地區的位置重要，175m水位蓄水以來，庫首區M≥2.0地震多發生于這個區域。綜合地層巖性、斷裂構造、滲透條件等，運用水-巖耦合數值模擬手段，探究該區域地震的成因機理是我們下一步研究的目標。