朱家正, 孫玉軍
中國(guó)地質(zhì)科學(xué)院地球深部探測(cè)中心, 北京 100037
川西地區(qū)位于我國(guó)南北地震帶的南段,是構(gòu)造形態(tài)復(fù)雜、地殼運(yùn)動(dòng)劇烈、地震活動(dòng)頻繁的地區(qū).在深大斷裂的控制下,塊體之間的水平運(yùn)動(dòng)和垂直差異運(yùn)動(dòng)都表現(xiàn)得十分明顯,而強(qiáng)烈地震大多數(shù)就發(fā)生在這些深大斷裂帶上(李碧雄等, 2014).同時(shí),我國(guó)西南地區(qū)水力資源豐富,近些年在該地區(qū)修建了多個(gè)大型水電站,如向家壩、溪洛渡、白鶴灘、烏東德等大型水電站,大型水電工程建設(shè)所面臨的水庫(kù)地震問(wèn)題是重要的科學(xué)問(wèn)題.一方面我國(guó)西南地區(qū)是構(gòu)造活動(dòng)十分強(qiáng)烈的區(qū)域,構(gòu)造地震活躍;另一方面,修建水庫(kù)后,大型水電工程蓄水造成的附加載荷和滲流作用一般都會(huì)改變地區(qū)原有的應(yīng)力狀態(tài).水庫(kù)地震一般震級(jí)較小,但是其造成的附加應(yīng)力改變可能與該區(qū)域初始的構(gòu)造應(yīng)力狀態(tài)疊加,從而造成大的構(gòu)造地震發(fā)生.如2008年發(fā)生在龍門山斷裂帶上的汶川地震,毫無(wú)疑問(wèn),如此大震級(jí)的地震發(fā)生是構(gòu)造運(yùn)動(dòng)的直接結(jié)果(Chen, 2009),但是后續(xù)的研究表明建立在龍門山斷裂帶上的紫坪鋪水庫(kù)對(duì)斷層的庫(kù)侖應(yīng)力影響為正值(Klose, 2012;周斌等, 2010;孫玉軍等, 2012;雷興林等, 2008;張貝和石耀霖, 2010),即對(duì)斷層滑動(dòng)有一定的促發(fā)作用,可能促使斷層滑動(dòng)提前發(fā)生(Ge et al., 2009).水庫(kù)蓄水可以誘發(fā)地震是一個(gè)已經(jīng)被證實(shí)的科學(xué)問(wèn)題(Gough and Gough, 1970;Gupta, 2002;Gupta and Chadha, 1995;Talwani, 1997;夏其發(fā), 1993;孫玉軍等, 2012).我國(guó)最著名的水庫(kù)地震是1962年發(fā)生在廣東的新豐江水庫(kù)地震,震級(jí)達(dá)到了6.1級(jí)(丁原章等, 1983;程惠紅等, 2012).因此在構(gòu)造活動(dòng)強(qiáng)烈的區(qū)域修建大型水庫(kù)等水電建設(shè)工程,其對(duì)該區(qū)域構(gòu)造活動(dòng)特別是大型斷裂的影響,是未來(lái)值得深入考慮和探究的科學(xué)問(wèn)題.該問(wèn)題不僅涉及到大型水電工程的安全穩(wěn)定性問(wèn)題,同時(shí)對(duì)建設(shè)區(qū)域的防震減災(zāi)也有重要的參考價(jià)值.
四川大崗山水電站修建于2013年,屬于大渡河上游修建的一系列水電站之一.由于龍門山斷裂、鮮水河斷裂、安寧河斷裂等一系列大型斷裂在此區(qū)域交匯,地質(zhì)背景復(fù)雜,構(gòu)造運(yùn)動(dòng)活躍,屬于地震高烈度區(qū).同時(shí),大崗山水庫(kù)屬于高山峽谷型水庫(kù),壩高達(dá)210 m,對(duì)于這類特殊性水庫(kù),其蓄水造成的附加災(zāi)害對(duì)該區(qū)重要的大型斷裂影響需要重點(diǎn)關(guān)注(程萬(wàn)正, 2013).
大崗山水庫(kù)庫(kù)區(qū)位于青藏高原東南緣川滇南北向構(gòu)造帶北端,川滇地塊、華南地塊、巴顏喀拉地塊的邊緣交界處(圖1),是地殼運(yùn)動(dòng)劇烈、構(gòu)造形態(tài)復(fù)雜、地震活動(dòng)十分頻繁的地區(qū)(張培震等, 2013;趙靜等, 2018).龍門山斷裂帶、鮮水河斷裂帶以及安寧河斷裂帶在此區(qū)域組成Y型構(gòu)造帶.大涼山斷裂、小金河斷裂、玉農(nóng)希斷裂等主要斷裂也在這一區(qū)域附近交匯,其中鮮水河斷裂南段的磨西斷裂以及大渡河斷裂都直接通過(guò)水庫(kù)蓄水區(qū)域(阮祥等, 2017).龍門山斷裂帶是一條長(zhǎng)約500 km、寬約30~50 km、沿NE-SW方向展布的巨大斷裂帶,按照由西向東的順序,龍門山斷裂帶主要包含龍門山后山斷裂(茂縣—汶川斷裂)、中央斷裂(映秀—北川斷裂)和山前斷裂(安縣—灌縣斷裂,亦稱彭縣—灌縣斷裂、江油—灌縣斷裂).這些斷裂(北西段)都以逆沖滑動(dòng)為主,兼具一定的右旋走滑分量(Xu et al., 2008;鄧起東等, 1994;陳國(guó)光等, 2007).從地震剖面上看,逆沖面呈鏟狀幾何形態(tài),上地殼向西陡傾,中地殼向西逐漸傾斜(Xu et al., 2008 ).張培震等的研究結(jié)果顯示,龍門山斷裂西南段存在(4±2.5)mm·a-1的縮短速率,以及(7.5±2)mm·a-1左旋走滑速率(張培震等, 2003).趙靜等采用DEFNODE負(fù)位錯(cuò)反演程序估算了龍門山斷裂西南端的閉鎖程度與變形狀態(tài),并綜合GPS反演結(jié)果和跨斷層水準(zhǔn)結(jié)果,認(rèn)為目前龍門山斷裂帶西南端在大部分段落處于強(qiáng)閉鎖狀態(tài),且依然有發(fā)生大地震的可能性(趙靜等, 2018).鮮水河斷裂帶北起甘孜東谷附近,向南經(jīng)過(guò)爐霍、道孚、康定一線,至石棉縣安順場(chǎng)一帶逐漸減弱消失(郭長(zhǎng)寶等, 2015).晚新生代以來(lái),鮮水河斷裂表現(xiàn)出強(qiáng)烈的左旋走滑運(yùn)動(dòng),是松潘—甘孜造山帶內(nèi)部一條大型走滑斷裂,橫切了松潘—甘孜造山帶的主體,系造山運(yùn)動(dòng)后期陸內(nèi)變形的產(chǎn)物,晚新生代以來(lái)的位移總規(guī)模在60 km左右(許志琴等, 2007).以乾寧惠遠(yuǎn)寺一帶為界,鮮水河斷裂帶北西段和南東段的活動(dòng)速率有所差異,北西段活動(dòng)速率為10~15 mm·a-1(聞學(xué)澤等, 1989),南東段小于10 mm·a-1,一般為 5 mm·a-1左右(郭長(zhǎng)寶等, 2015).李鐵明等通過(guò)重力、GPS資料計(jì)算結(jié)果,得出鮮水河斷裂磨西段先進(jìn)整體呈現(xiàn)左旋走滑,速率為4.41 mm·a-1(李鐵明等, 2019).
圖1 大崗山水庫(kù)區(qū)域構(gòu)造背景(a)及庫(kù)區(qū)地形圖(b)(a) 中橙色實(shí)線代表二級(jí)地塊邊界,藍(lán)色實(shí)線代表一級(jí)地塊邊界,黑色實(shí)線代表主要斷裂,黃色方框是研究區(qū)域所在位置,黃色五角星代表水庫(kù)大壩位置; (b) 中藍(lán)色填充區(qū)域?yàn)樗畮?kù)主要蓄水區(qū)域,黃色方塊代表大壩所在位置.Fig.1 Regional tectonic background (a) and topographic map (b) of Dagangshan ReservoirThe solid orange lines in the figure (a) represent the boundaries of the secondary blocks. The solid blue lines represent the boundaries of the primary blocks. The solid black lines represent the main faults. The yellow box represents the location of the study area, and the yellow pentacle represents the location of the reservoir dam. The blue filled area in the figure (b) is the main water storage area of the reservoir, and the yellow square represents the location of the dam.
大崗山水庫(kù)在2013年建成后,大致分兩階段使得蓄水水位達(dá)到預(yù)定水位.于2014年3月開(kāi)啟第一階段蓄水,水位達(dá)到1000 m;于2014年11月開(kāi)啟第二階段蓄水,于2015年5月達(dá)到預(yù)定水位,水庫(kù)水位由2013年蓄水前的960 m水位上升到1130 m的正常蓄水位.庫(kù)區(qū)的地震數(shù)據(jù)表明(圖2),水庫(kù)蓄水后研究區(qū)域內(nèi)地震發(fā)生頻率隨著水位的變化發(fā)生了明顯的增加,特別是伴隨著水位的兩次大幅度變化,庫(kù)區(qū)周邊地震頻度以及震級(jí)均發(fā)生增大.本研究將101.9°E—102.4°E, 29.2°N—29.8°N作為研究范圍,以國(guó)家地震臺(tái)網(wǎng)的地震數(shù)據(jù)為依據(jù),繪出蓄水前后三年內(nèi)庫(kù)區(qū)地震平面分布圖(圖3).可以看到,在水庫(kù)蓄水之前,地震主要分布在庫(kù)區(qū)西側(cè)鮮水河斷裂,多為ML3級(jí)以下微震,在鮮水河斷裂西側(cè)呈條帶狀分布;2015年5月水庫(kù)經(jīng)過(guò)第二階段蓄水達(dá)到正常水位后,該區(qū)域的地震分布仍然主要分布在庫(kù)區(qū)西北的鮮水河斷裂磨西段,但叢集性更強(qiáng),且出現(xiàn)ML3~4級(jí)以上地震8次,ML4~5級(jí)地震3次,活動(dòng)性變強(qiáng).但在庫(kù)區(qū)東部的龍門山斷裂南段,在蓄水達(dá)到正常水位后幾乎不再出現(xiàn)微震,庫(kù)區(qū)南側(cè)的大涼山斷裂附近的地震活動(dòng)性也有一定的降低.該現(xiàn)象是否與水庫(kù)蓄水有直接關(guān)系值得深入研究.
圖2 研究區(qū)域地震活動(dòng)性與水庫(kù)水位關(guān)系Fig.2 The relationship between seismic activity and reservoir water level in the study area
圖3 水庫(kù)蓄水前后研究區(qū)域地震活動(dòng)情況(a) 2011年11月—2014年11月庫(kù)區(qū)地震分布; (b) 2014年11月—2017年11月庫(kù)區(qū)地震分布. 黃色五角星標(biāo)識(shí)為水庫(kù)所在位置.a, b, c, d分別代表龍門山斷裂、鮮水河斷裂、大涼山斷裂以及錦屏山—小金河斷裂.○代表震中位置,綠色為ML0~2級(jí),黃色為ML2~3級(jí),橙色代表ML3~4級(jí),紅色代表ML4級(jí)以上地震.紅色橢圓所圈出區(qū)域?yàn)樾钏蟮卣鸹顒?dòng)性增強(qiáng)區(qū)域,藍(lán)色橢圓所圈出的區(qū)域?yàn)樾钏蟮卣鸹顒?dòng)性降低區(qū)域.Fig.3 Seismic activity in the study area before and after the reservoir storage(a) The distribution of earthquakes in the reservoir area from November 2011 to November 2014; (b) The distribution of earthquakes in the reservoir area from November 2014 to November 2017. The yellow star indicates the location of the dam. a, b, c, d represent Longmenshan fault, Xianshuihe fault, Daliangshan fault and Jinpingshan-Xiaojinhe fault respectively. ○ represents the location of the epicenter. Green represents magnitude ML0~2, yellow represents magnitude ML2~3, orange represents magnitude ML3~4, and red represents above magnitude ML4. The red ellipse circles the area of increased seismicity after water storage, and the blue ellipse circles the area of decreased seismicity after water storage.
為了定量分析大崗山水庫(kù)蓄水對(duì)該區(qū)域斷層的影響,本文依據(jù)構(gòu)造地質(zhì)特征及精細(xì)的DEM數(shù)據(jù)建立了三維孔隙彈性有限元數(shù)值模型(圖4),通過(guò)計(jì)算大崗山水庫(kù)蓄水過(guò)程造成庫(kù)區(qū)斷層的庫(kù)侖應(yīng)力變化,定量分析了大崗山水庫(kù)蓄水對(duì)地震活動(dòng)性的影響.
圖4 三維孔隙彈性有限元模型(模型中不同顏色代表不同材料參數(shù))Fig.4 Three-dimensional poroelastic finite element model (Different colors in the model represent different material parameters)
孔隙彈性體的一般本構(gòu)關(guān)系為
(1)
水的滲流擴(kuò)散方程(暫不考慮飽和孔隙水不排水效應(yīng)的影響)為(Bell and Nur,1978;Roeloffs, 1988)
(2)
(3)
其中c為水力學(xué)擴(kuò)散系數(shù),單位m2·s-1,k是介質(zhì)的滲透率,單位m2,η是水的流體動(dòng)力學(xué)黏滯性系數(shù),單位Pa·s,泊松比、剪切模量、Skempton系數(shù)對(duì)結(jié)果的影響主要是通過(guò)改變擴(kuò)散系數(shù),本文將c作為計(jì)算參數(shù).c的取值范圍主要參考Talwani等(2007).
在小形變的幾何方程和平衡方程的約束下,施加合理的邊界條件就可以進(jìn)行有限元求解.
幾何方程:
(4)
平衡方程(忽略體力項(xiàng)):
(5)
本文采用孫玉軍等開(kāi)發(fā)的三維孔隙彈性有限元模擬程序進(jìn)行計(jì)算(孫玉軍等, 2012).流體滲流和孔隙彈性耦合的研究方法不僅可以用于研究水庫(kù)地震,同時(shí)在研究頁(yè)巖氣開(kāi)采過(guò)程中壓裂誘發(fā)地震時(shí)也被廣泛應(yīng)用(Hemam et.al., 2021; Wong et.al., 2021; Yang et.al., 2021).
三維模型的建立考慮地層的垂向分層,模型表面層采用實(shí)際的地表地形數(shù)據(jù),本研究采用SRTM30 m數(shù)字高程數(shù)據(jù)(https:∥www.usgs.gov/products/maps/topo-maps).設(shè)定模型的深度為25 km,模型范圍為101.9°E—102.4°E,29.2°N—29.8°N,采用三棱柱單元進(jìn)行網(wǎng)格剖分,并對(duì)庫(kù)區(qū)部分進(jìn)行網(wǎng)格加密,模型總結(jié)點(diǎn)數(shù)702259,總單元數(shù)1399412.數(shù)值模型可以考慮實(shí)際的分層性,但是由于缺乏實(shí)際滲流參數(shù),本研究只區(qū)分了斷層和非斷層區(qū)域,對(duì)其設(shè)置了不同的模型參數(shù).水平層狀網(wǎng)格中,對(duì)蓄水區(qū)域進(jìn)行加密,分辨率為60 m,非蓄水區(qū)域網(wǎng)格大小為900 m.垂向上網(wǎng)格分為15層,深度為25 km,底部四層網(wǎng)格無(wú)斷層穿過(guò),斷層產(chǎn)狀傾向北西西,傾角83°.
在研究區(qū)域范圍內(nèi),對(duì)水庫(kù)蓄水位以上的其他區(qū)域,給定力的邊界和孔隙壓力邊界均為0,對(duì)水位以下地區(qū)給定力的邊界和孔隙壓力邊界條件為
p(x,y,z,t)=ρωgh(x,y,z,t),
(6)
ρω為水的密度,g=9.8 m·s-2為重力加速度,h(x,y,z,t)為t時(shí)刻水庫(kù)水位減去不同地點(diǎn)的地形高程值,即實(shí)際水深,p(x,y,z,t)為對(duì)應(yīng)的水庫(kù)蓄水造成的庫(kù)底水壓力,對(duì)于力的邊界,模型側(cè)邊界和底邊界均為滑動(dòng)邊界條件,對(duì)于孔隙壓力,側(cè)面給定為等滲流梯度,底面固定孔隙壓力為0.
庫(kù)侖應(yīng)力的變化量定義為(Harris,1998;石耀霖和曹建玲, 2010)
Δcfs=Δτ+μ(Δσn+Δp),
(7)
Δτ為所考慮的斷層面上剪應(yīng)力的變化量,μ是摩擦系數(shù),Δσn為斷層面上的正應(yīng)力改變量(拉伸為正),Δp為孔隙壓力改變量.本文以coulomb-s代表庫(kù)侖應(yīng)力變化量.
同時(shí),為了比較不同參數(shù)對(duì)計(jì)算結(jié)果的影響,分別建立不同的模型進(jìn)行比較(表1).
表1 各個(gè)模型中計(jì)算參數(shù)的選取Table 1 Selection of calculation parameters in each model
鮮水河斷裂南東段的磨西斷裂主要走向?yàn)?60°(李大虎等, 2015),在靠近庫(kù)區(qū)的位置取其走向?yàn)?66°,根據(jù)程佳、徐錫偉等對(duì)鮮水河南端水平和垂直滑動(dòng)速率的總結(jié),將斷層滑動(dòng)角參數(shù)設(shè)置為21°(徐錫偉等, 2003;程佳等, 2012).龍門山斷裂南段大渡河斷裂的斷層參數(shù)參考前人的區(qū)域應(yīng)力場(chǎng)及震源機(jī)制解研究來(lái)選取(易桂喜等, 2013;李君等, 2019;杜瑤等, 2016;聞學(xué)澤等, 2008).計(jì)算庫(kù)侖應(yīng)力所采用的接收斷層參數(shù)選取如表2所示.
表2 計(jì)算庫(kù)侖應(yīng)力采用的各斷層參數(shù)Table 2 The parameters of each fault used to calculate the change of coulomb stress (coulomb-s)
為討論不同模型參數(shù)對(duì)計(jì)算結(jié)果影響,在水庫(kù)蓄水后庫(kù)區(qū)西北側(cè)選取參考點(diǎn)1坐標(biāo)為(102.12°,29.53°,-5 km),分別對(duì)比三組模型在該點(diǎn)處的孔隙壓力以及庫(kù)侖應(yīng)力變化,庫(kù)侖應(yīng)力的計(jì)算采用鮮水河斷裂磨西段的斷層參數(shù)(表2).圖5中coulomb-s代表庫(kù)侖應(yīng)力值.
圖5a、5b為模型1計(jì)算得到結(jié)果.模型1中通過(guò)保持其他參數(shù)不變而只改變了斷層的彈性模量以探究其對(duì)結(jié)果影響.從圖中可以看出,在2013年12月水位第一次上漲時(shí)水庫(kù)蓄水就對(duì)該參考點(diǎn)位置造成正的庫(kù)侖應(yīng)力影響,但該蓄水過(guò)程中孔隙壓力幾乎沒(méi)有變化,說(shuō)明此時(shí)水庫(kù)蓄水的彈性載荷效應(yīng)在起主要作用,隨著水庫(kù)在2014年12月水位快速增長(zhǎng)110 m,三種不同模型參數(shù)對(duì)應(yīng)的庫(kù)侖應(yīng)力變化量也增加明顯,同時(shí)參考點(diǎn)的庫(kù)侖應(yīng)力隨斷層彈性模量的增大而減小,但由于三個(gè)模型并未改變滲流狀態(tài),因此三種模型所對(duì)應(yīng)的孔隙壓力基本一樣(圖5b),隨著后期孔隙壓力的增加,水庫(kù)蓄水過(guò)程造成了正的庫(kù)侖應(yīng)力增加.
圖5c、5d為模型2計(jì)算結(jié)果.模型2只改變了斷層的擴(kuò)散系數(shù)而保持其他系數(shù)不變.可以從結(jié)果中看出孔隙壓力的結(jié)果在2015年4月出現(xiàn)明顯改變;而庫(kù)侖應(yīng)力的三條結(jié)果曲線在2014年12月水位大幅度改變之前基本保持一致,計(jì)算結(jié)果顯示此時(shí)水的滲流還未達(dá)到該深度.在2014年12月之后庫(kù)侖應(yīng)力隨著孔隙壓力的增大而逐漸變大,同時(shí)斷層的擴(kuò)散系數(shù)越大,在參考點(diǎn)處造成的庫(kù)侖應(yīng)力變化量則越小.
圖5e、5f為保持?jǐn)鄬拥臄U(kuò)散系數(shù)不變,而改變非斷層區(qū)域的擴(kuò)散系數(shù).結(jié)果顯示非斷層區(qū)域的擴(kuò)散系數(shù)對(duì)結(jié)果影響較為明顯,在第二階段蓄水之后,當(dāng)擴(kuò)散系數(shù)c=1 m2·s-1時(shí),參考點(diǎn)處的庫(kù)侖應(yīng)力和孔隙壓力值都快速增長(zhǎng),而擴(kuò)散系數(shù)低于該值的另外兩個(gè)模型顯示,庫(kù)侖應(yīng)力和孔隙壓力變化較為平緩.
圖5g、5h為改變非斷層區(qū)域的排水泊松比所得到的計(jì)算結(jié)果.計(jì)算結(jié)果顯示,在保持?jǐn)U散系數(shù)不變的情況下,改變非斷層區(qū)域的排水泊松比對(duì)孔隙壓力的結(jié)果沒(méi)有影響.水庫(kù)蓄水初期庫(kù)侖應(yīng)力隨著排水泊松比的增大而減小,隨著時(shí)間的增加三種計(jì)算結(jié)果在一點(diǎn)交會(huì)后,兩者關(guān)系轉(zhuǎn)為正相關(guān).該結(jié)果可以通過(guò)孔隙彈性的本構(gòu)關(guān)系來(lái)理解,也即單純的彈性介質(zhì)下,庫(kù)侖應(yīng)力隨排水泊松比增大而減小,但是根據(jù)孔隙彈性的本構(gòu)關(guān)系,孔隙壓力的增加,高排水泊松比更有利于庫(kù)侖應(yīng)力的增大.
圖5 不同模型計(jì)算得到的庫(kù)侖應(yīng)力變化量以及孔隙壓力隨時(shí)間變化(a、b) 模型1結(jié)果; (c、d) 模型2結(jié)果; (e、f) 模型3結(jié)果. (g、h) 模型4結(jié)果; (b、h)中紅色、綠色和粉色線重合.Fig.5 The change of Coulomb stress and pore pressure calculated by different models vary with time(a,b) The results of Model 1; (c,d) The results of Model 2; (e,f) The results of Model 3. The red, green and pink lines in (b) coincide.
以上幾種模型中都可以看到在水庫(kù)水位變化引起庫(kù)侖應(yīng)力變化主要受到兩種效應(yīng)影響,一種是蓄水造成彈性載荷變化引起的庫(kù)侖應(yīng)力變化,另一種是水的滲流造成孔隙壓力變化而引起的庫(kù)侖應(yīng)力變化.但這兩種效應(yīng)在時(shí)間上的表現(xiàn)有所不同,蓄水造成彈性載荷變化所引起的庫(kù)侖應(yīng)力變化是瞬時(shí)的,如第一階段蓄水之后庫(kù)侖應(yīng)力都隨之快速增加,此時(shí)斷層彈性模量對(duì)計(jì)算結(jié)果有較大影響,但此時(shí)孔隙壓力并未明顯改變.而第二種效應(yīng)中孔隙壓力改變所造成的庫(kù)侖應(yīng)力變化相對(duì)于水庫(kù)蓄水過(guò)程有所滯后,如第二階段蓄水之后的庫(kù)侖應(yīng)力增長(zhǎng)都隨著孔隙壓力的變化而變化,說(shuō)明水庫(kù)達(dá)到設(shè)計(jì)水位后,主要是水的滲流作用影響庫(kù)侖應(yīng)力,而對(duì)此影響較大的參數(shù)是斷層或非斷層區(qū)域的擴(kuò)散系數(shù),擴(kuò)散系數(shù)越大,孔隙壓力越大,從而對(duì)參考點(diǎn)造成的庫(kù)侖應(yīng)力改變?cè)酱?
從孔隙彈性介質(zhì)材料來(lái)看,水庫(kù)蓄水過(guò)程所產(chǎn)生的這兩種效應(yīng),可能會(huì)互相疊加增強(qiáng),也有可能會(huì)互相抵消減弱.為了更詳細(xì)的分析兩種效應(yīng)的作用,我們?cè)谟?jì)算庫(kù)侖應(yīng)力的時(shí)候可以分別對(duì)其進(jìn)行考慮,即考慮水庫(kù)滲流造成的庫(kù)侖應(yīng)力改變(彈性和滲流的綜合效應(yīng))和不考慮水庫(kù)滲流造成的庫(kù)侖應(yīng)力變化(彈性效應(yīng)).
綜合上述三種模型,以及前人的研究成果(Talwani et al., 2007;陳建業(yè)等, 2011)選取模型5為最終的模型計(jì)算參數(shù).在庫(kù)區(qū)西側(cè)的磨西斷裂選取參考點(diǎn)1(102.12°E,29.53°N,-5 km),大渡河斷裂選取參考點(diǎn)2(102.20°E,29.57°N,-5 km)(圖3a),對(duì)水庫(kù)蓄水后的地震活動(dòng)性增強(qiáng)和減弱區(qū)域進(jìn)行庫(kù)侖應(yīng)力和孔隙壓力進(jìn)行分析.對(duì)點(diǎn)1庫(kù)侖應(yīng)力計(jì)算選取磨西斷裂的斷層參數(shù),對(duì)點(diǎn)2的庫(kù)侖應(yīng)力計(jì)算采用大渡河斷裂的斷層參數(shù)(表2).圖中coulomb-e為不考慮孔隙壓力影響計(jì)算得到的庫(kù)侖應(yīng)力變化量.
對(duì)比兩圖發(fā)現(xiàn),圖6a中在水庫(kù)蓄水后不含孔隙壓力的庫(kù)侖應(yīng)力coulomb-e(彈性效應(yīng))為正值,并且隨水位增加而迅速增長(zhǎng),使斷層更加危險(xiǎn);而圖6b的結(jié)果與之相反,不含孔隙壓力的庫(kù)侖應(yīng)力coulomb-e(彈性效應(yīng))為負(fù)值,隨水位增長(zhǎng)而快速下降,對(duì)斷層起到抑制作用.這是由于水庫(kù)蓄水區(qū)域跨越了龍門山和鮮水河兩大斷裂的南端,位置分別處于磨西斷裂的下盤,以及大渡河斷裂的上盤.因此庫(kù)水對(duì)斷層施加的彈性載荷使得磨西斷裂的主要區(qū)域地震活動(dòng)性增加,使得大渡河斷裂的主要控制區(qū)域地震活動(dòng)性減弱.但是,隨著水位保持高位,水不斷向斷層中滲流,孔隙壓力隨之增大,使得庫(kù)區(qū)兩側(cè)的庫(kù)侖應(yīng)力coulomb-s都轉(zhuǎn)變?yōu)檎?斷層趨向于更加危險(xiǎn).
圖6 地震活動(dòng)性增強(qiáng)區(qū)域與減弱區(qū)域選定點(diǎn)的庫(kù)侖應(yīng)力和孔隙壓力變化(a) 地震活動(dòng)性增強(qiáng)區(qū)域點(diǎn)1結(jié)果; (b) 地震活動(dòng)性減弱區(qū)域點(diǎn)2結(jié)果.Fig.6 Changes in Coulomb stress and pore pressure at selected points in the areas of enhanced and weakened seismic activity(a) The result of point 1 in the seismic activity enhancement area; (b) The result of point 2 in the seismic activity weakened area.
將磨西斷裂參數(shù)作為計(jì)算庫(kù)侖應(yīng)力的參數(shù),得到2017年1月的庫(kù)侖應(yīng)力變化量與孔隙壓力的平面分布(圖7).孔隙壓力增大的區(qū)域?qū)?yīng)水庫(kù)主蓄水區(qū)域,且隨著深度的增加,孔隙壓力值逐漸減小,如在深度5 km孔隙壓力最大值為63 kPa,而在10 km和15 km孔隙壓力最大值分別為22 kPa和10 kPa;同時(shí)由于斷層的擴(kuò)散系數(shù)一般較大,孔隙壓力在斷層區(qū)域一般較高.庫(kù)侖應(yīng)力變化量正值區(qū)域除了在水庫(kù)庫(kù)區(qū)范圍外,主要分布在庫(kù)區(qū)的北西和南東兩個(gè)方向,同時(shí)隨深度增加逐漸減小,在5 km深度處庫(kù)侖應(yīng)力變化量正極大值可達(dá)45 kPa,在10 km和15 km庫(kù)侖應(yīng)力變化量正極大值分別為17 kPa和7 kPa.庫(kù)侖應(yīng)力變化量的負(fù)值區(qū)域主要分布在庫(kù)區(qū)東北側(cè)和西南側(cè),在5 km深度負(fù)的極大值可達(dá)-6 kPa,在10 km和15 km負(fù)的極大值分別為-4 kPa和-2 kPa.深度為10 km以上時(shí),水的滲流造成的孔隙壓力較小,庫(kù)侖應(yīng)力變化量的區(qū)域性分布更為明顯,庫(kù)區(qū)西南方向出現(xiàn)了庫(kù)侖應(yīng)力變化量的負(fù)值區(qū)域,在水庫(kù)蓄水后的地震增加區(qū)域(庫(kù)區(qū)北西側(cè))看到明顯的2 kPa左右的正值區(qū)域,說(shuō)明水庫(kù)蓄水對(duì)該區(qū)域地震的發(fā)生具有促進(jìn)作用.
圖7 2017年1月模型計(jì)算得到的深度為5 km (a, b),10 km (c, d),15 km (e, f)的庫(kù)侖應(yīng)力變化量以及孔隙壓力分布(單位:Pa)○ 代表第二階段蓄水之后距離剖面垂直距離范圍2 km內(nèi)的地震震中位置在平面上的投影,灰色代表震級(jí)ML2以下微震,黃色代表ML2~3級(jí)地震,橙色代表ML3~4級(jí)地震,紅色代表ML4級(jí)以上地震.Fig.7 The coulomb stress and pore pressure distribution at the depth of 5 km (a, b),10 km (c, d), and 15 km (e, f) in January, 2017 (unit: Pa)○ Represents the projection of the earthquake epicenter position on the plane within 2 kilometers of the vertical distance from the profile after the second stage of water storage. Gray represents seismic magnitude below ML2, yellow represents magnitude ML2~3, orange represents magnitude ML3~4, and red represents magnitude above ML4.
以參考點(diǎn)(102.12°,29.53°)為中心位置對(duì)計(jì)算結(jié)果做縱向剖面(圖8),孔隙壓力主要分布在水庫(kù)區(qū)域且有隨著斷層在地下的延展更快速地?cái)U(kuò)散(圖8c),庫(kù)侖應(yīng)力變化量的分布與圖7對(duì)應(yīng),在庫(kù)區(qū)南北兩側(cè)形成負(fù)值區(qū)域,中心位置為正值;庫(kù)區(qū)東側(cè)以及南側(cè)的5~10 km深處存在庫(kù)侖應(yīng)力變化量的負(fù)值區(qū)域(圖8b、8d),西側(cè)沿?cái)鄬拥难诱狗植紟?kù)侖應(yīng)力變化量正值區(qū)域.
圖8 模型4計(jì)算得到的2017年1月不同剖面孔隙壓力(a,c)和庫(kù)侖應(yīng)力變化量(b, d)分布圖剖面A-A′為102.15°E的縱切面,剖面B-B′為29.53°N的橫切面. 剖面位置已在圖3a中顯示.Fig.8 The pore pressure (a, c) and Coulomb stress (b, d) distribution of different profiles calculated by Model 4 in January 2017 (The positions of the two profiles are shown in Fig.3a)
從應(yīng)力觸發(fā)區(qū)域和地震活動(dòng)的對(duì)應(yīng)關(guān)系上,可以看到水庫(kù)蓄水對(duì)10 km和15 km深度范圍的地震在空間位置上具有明顯的對(duì)應(yīng)關(guān)系:在10 km深度處的應(yīng)力影響區(qū)邊緣出現(xiàn)了較為集中的ML2~3級(jí)地震分布,15 km深度范圍內(nèi)的ML3~5級(jí)的地震也集中分布在西北側(cè)應(yīng)力觸發(fā)區(qū)域的邊緣,說(shuō)明水庫(kù)蓄水對(duì)該區(qū)域地震發(fā)生具有促進(jìn)作用,而不同深部的應(yīng)力影區(qū)域的地震活動(dòng)均受到抑制.
以大渡河斷裂參數(shù)(表2)作為計(jì)算依據(jù)得出2017年1月庫(kù)侖應(yīng)力變化量(圖9a, 9c)及不考慮水滲流的庫(kù)侖應(yīng)力變化量(圖9b, 9d)平面分布圖.不考慮水滲流的庫(kù)侖應(yīng)力變化量在庫(kù)區(qū)東側(cè)存在一個(gè)明顯的負(fù)值區(qū)域,在5 km深度靠近水庫(kù)大壩位置的庫(kù)侖應(yīng)力變化量負(fù)極大值達(dá)到-12 kPa,考慮水滲流的情況下該值增加至-8 kPa左右,且與采用磨西斷裂參數(shù)時(shí)一樣,在蓄水區(qū)域出現(xiàn)正值區(qū),5 km深度處庫(kù)侖應(yīng)力變化量正值極大值達(dá)到50 kPa,10 km深度處的庫(kù)侖應(yīng)力變化量分布與采用磨西斷裂得到的結(jié)果相比,庫(kù)侖應(yīng)力變化量(考慮水滲流)負(fù)的極值增大1 kPa左右,達(dá)到-5.6 kPa.對(duì)比采用磨西斷裂參數(shù)計(jì)算得到的結(jié)果可以看出,盡管采用的兩條斷層參數(shù)有所改變,但水庫(kù)蓄水所造成的庫(kù)侖應(yīng)力變化量正、負(fù)區(qū)域(即應(yīng)力觸發(fā)和應(yīng)力影區(qū)域)并未明顯改變,除了造成庫(kù)區(qū)的庫(kù)侖應(yīng)力變化量為正值外,在庫(kù)區(qū)北西和南東為庫(kù)侖應(yīng)力變化量正值區(qū)(應(yīng)力觸發(fā)區(qū)域),庫(kù)區(qū)北西側(cè)相對(duì)較大;而在庫(kù)區(qū)北東和南西則為庫(kù)侖應(yīng)力變化量負(fù)值區(qū)(應(yīng)力影區(qū)域),庫(kù)區(qū)北東負(fù)值更明顯.該影響與水庫(kù)蓄水后地震活動(dòng)性比較一致.
圖9 模型4計(jì)算得到的2017年1月深度為5 km(a, b)和10 km (c, d)的庫(kù)侖應(yīng)力變化量分布 (a, c)和(b, d)分別為考慮水滲流和不考慮水滲流得到了庫(kù)侖應(yīng)力變化量分布(單位:Pa). ○代表第二階段蓄水之后距離剖面垂直距離范圍2 km內(nèi)的地震震中位置在平面上的投影,灰色為ML2級(jí)以下,黃色為ML2~3級(jí),橙色代表ML3~4級(jí),紅色代表ML4級(jí)以上地震.藍(lán)色三角形代表石棉地震震中投影.Fig.9 The distribution of the change of coulomb stress at the depth of 5 km (a, b) and 10 km (c, d) in January, 2017 for the model 4 (unit: Pa)(a, c) and (b, d) are the change of Coulomb stress with and without considering pore-pressure, respectively. ○ Represents the projection of the earthquake epicenter position on the plane within 2 kilometers of the vertical distance from the profile after the second stage of water storage. Gray represents microseisms below magnitude ML2, yellow represents magnitude ML2~3 earthquakes, orange represents magnitude ML3~4 earthquakes, and red represents earthquakes above magnitude ML4. The blue triangle represents the asbestos earthquake epicenter projection.
本文采用三維孔隙彈性模型,研究大崗山水庫(kù)蓄水前后水位變化對(duì)庫(kù)區(qū)地震及斷層的影響.根據(jù)地震臺(tái)網(wǎng)數(shù)據(jù),已知經(jīng)過(guò)水位變化的第二階段蓄水之后,水庫(kù)附近的微震活動(dòng)變化明顯,首先,庫(kù)區(qū)西北的鮮水河斷裂南段附近地震震級(jí)增大,出現(xiàn)了ML3級(jí)以上的地震,且數(shù)量增多,叢集性更明顯,由于蓄水位置處于磨西斷裂下盤,促進(jìn)了斷層活動(dòng)性,根據(jù)計(jì)算結(jié)果,在地震叢集增強(qiáng)的區(qū)域靠近庫(kù)區(qū)的磨西斷裂一側(cè)為庫(kù)侖應(yīng)力變化量正值區(qū)(應(yīng)力觸發(fā)區(qū));而在地震活動(dòng)性降低的庫(kù)區(qū)北東區(qū)域,庫(kù)侖應(yīng)力變化量為負(fù)值(應(yīng)力影區(qū)).但是隨著滲流的持續(xù),庫(kù)侖應(yīng)力變化量負(fù)值區(qū)域也會(huì)逐漸變?yōu)檎?地震活動(dòng)性增強(qiáng).根據(jù)前人研究成果,鮮水河斷裂帶的最大剪切應(yīng)變速率為40~60×10-9/a(Wang and Shen, 2020).如果取最大剪切應(yīng)變速率為50×10-9/a,假設(shè)彈性模量為80 GPa的情況下,應(yīng)力累積速率為4 kPa/a,而蓄水造成的庫(kù)侖應(yīng)力變化量增加在淺層5 km深度最大可達(dá)到45 kPa(圖10),這相當(dāng)于11年的構(gòu)造運(yùn)動(dòng)產(chǎn)生的應(yīng)力累積.
地震增強(qiáng)區(qū)域和減弱區(qū)域的平面方向與臨近的鮮水河斷裂和龍門山斷裂存在一定的夾角,這是否對(duì)應(yīng)了兩條斷裂帶地下的展布狀態(tài),需要進(jìn)一步討論.由于水庫(kù)位于磨西斷裂的下盤且主要蓄水區(qū)域都被其直接穿過(guò),導(dǎo)致無(wú)論水的滲流還是庫(kù)水的彈性載荷都促使其更加活躍;而大渡河斷裂的上盤由于受到水庫(kù)的垂向重力載荷,導(dǎo)致活動(dòng)性減弱,地震減少.但是隨著水的不斷滲流,以及水位的穩(wěn)定,水庫(kù)彈性載荷效應(yīng)造成的對(duì)地震的抑制作用會(huì)不斷被削弱,庫(kù)區(qū)的庫(kù)侖應(yīng)力可能都處于增加狀態(tài),即整個(gè)區(qū)域的地震活動(dòng)性更強(qiáng).根據(jù)地震數(shù)據(jù),該區(qū)域在蓄水之后三年的時(shí)間中都沒(méi)有再發(fā)生地震,說(shuō)明現(xiàn)實(shí)中庫(kù)水的彈性載荷對(duì)該區(qū)域地震的抑制效應(yīng)比預(yù)計(jì)的更強(qiáng);另外一種可能是水庫(kù)區(qū)域地殼結(jié)構(gòu)的黏彈性效應(yīng),黏彈性效應(yīng)往往也造成庫(kù)侖應(yīng)力分布的不均勻,也可能對(duì)該區(qū)域造成負(fù)的庫(kù)侖應(yīng)力,從而造成滲流效應(yīng)在該區(qū)域發(fā)生作用需要更長(zhǎng)的時(shí)間.因此,在未來(lái)的研究中對(duì)于地殼結(jié)構(gòu)的黏彈性效應(yīng)也應(yīng)該有所考慮.水庫(kù)北東側(cè)的斷裂活動(dòng)受到抑制,產(chǎn)生閉鎖效應(yīng),造成龍門山斷裂西南段的應(yīng)變累積速率加快,水位的快速下降有可能造成此前積累的應(yīng)變能快速釋放,增加大震發(fā)生的可能性.
由于模型和實(shí)際斷層詳細(xì)數(shù)據(jù)的限制,沒(méi)有給出鮮水河和龍門山斷裂的鏟狀結(jié)構(gòu),實(shí)際情況中斷層的傾角可能會(huì)隨著深度增加而減小,因此水的滲流影響會(huì)在深部隨著斷層擴(kuò)散至更遠(yuǎn)的區(qū)域,造成西北部區(qū)域的孔隙壓力增加值比模型計(jì)算得到的值更大,該區(qū)域的斷層可能更加活躍.另外,在鮮水河斷裂的西北區(qū)域也可能存在與斷裂平行的小斷裂,如果這些小斷裂與鮮水河主斷裂參數(shù)一致,那么目前計(jì)算得到的庫(kù)侖應(yīng)力也適用于這些小斷裂,從而促進(jìn)斷裂活動(dòng).因此未來(lái)如果獲取了該區(qū)域更為詳細(xì)的活動(dòng)斷裂探測(cè)數(shù)據(jù),則可以利用該模型開(kāi)展更為詳細(xì)的計(jì)算分析.
本文的主要結(jié)論如下:
(1)從孔隙彈性材料的角度看,水庫(kù)蓄水所造成的庫(kù)侖應(yīng)力改變主要有兩種效應(yīng),一種是蓄水造成彈性載荷變化引起的庫(kù)侖應(yīng)力變化,另一種是水的滲流造成孔隙壓力變化而引起的庫(kù)侖應(yīng)力變化.這兩種效應(yīng)在時(shí)間上的表現(xiàn)有所不同,蓄水造成彈性載荷變化所引起的庫(kù)侖應(yīng)力變化是瞬時(shí)的;而孔隙壓力改變所造成的庫(kù)侖應(yīng)力變化相對(duì)于水庫(kù)蓄水過(guò)程有所滯后.
(2)依據(jù)斷層實(shí)際滑動(dòng)參數(shù)計(jì)算得到的庫(kù)侖應(yīng)力結(jié)果來(lái)看,大崗山水庫(kù)蓄水過(guò)程所造成的庫(kù)侖應(yīng)力改變使得鮮水河斷裂磨西段庫(kù)侖應(yīng)力增加,而使得大渡河斷裂庫(kù)侖應(yīng)力降低.這與水庫(kù)蓄水后地震臺(tái)網(wǎng)觀測(cè)得到的磨西斷裂西側(cè)地震活動(dòng)性明顯增強(qiáng),而大渡河斷裂東側(cè)地震活動(dòng)性減弱的事實(shí)一致.
致謝感謝兩位匿名審稿人和編輯提出的建設(shè)性意見(jiàn).