胡幸平, 崔效鋒*, 張廣偉, 王甘嬌, Arno Zang,史丙新, 姜大偉
1 應急管理部國家自然災害防治研究院, 北京 100085 2 江西省地震局, 南昌 330039 3 GFZ German Research Centre for Geosciences, Telegrafenberg, 14473 Potsdam, Germany 4 四川省地震局, 成都 610041
四川長寧地區位于四川盆地南緣,揚子板塊西緣,區域構造包括四川盆地與蜀南低陡褶皺帶(王適擇等, 2013),是我國主要的井礦鹽和頁巖氣產區之一,存在多處采鹽、采氣及廢水回注井(阮祥等, 2008; 王玉滿等, 2016).長寧地區地震頻發且近年來增幅顯著,據中國地震臺網中心的數據,該地區在2018年12月以前并沒有5級以上地震,而從2018年12月16日興文MS5.7地震發生以來,陸續發生了多次5級以上地震,包括2019年1月3日珙縣MS5.1地震,以及2019年6月17日長寧MS6.0地震及其后20天內的4次5級以上中強震(圖1).由于毗鄰鹽氣開采區,且注水誘發地震又是國內外長期以來的研究熱點(Healy et al., 1968; 刁守中等, 1987; Segal, 1989; 張致偉等, 2012; Ellsworth, 2013; Bao and Eaton, 2016; Zang et al., 2019),因此,長寧地區的地震活動及其成因機理引起了科學界乃至整個社會的廣泛關注(阮祥等, 2008; 朱航和何暢, 2014; Sun et al., 2017; Lei et al., 2017, 2019; Meng et al., 2019;何登發等, 2019; 易桂喜等, 2019; 劉敬光等, 2019; 胡曉輝等, 2020; 常祖峰等, 2020).
針對長寧地區地震活動的已有研究指出,2019年6月17日長寧MS6.0地震及其余震序列發生于長寧背斜這一當地主要地質構造的軸部,其震源機制解的P軸主要沿北東東-南西西向,整體呈逆沖兼走滑的構造特征(易桂喜等, 2019);而2018年12月16日興文MS5.7地震和2019年1月3日珙縣MS5.1地震發生于長寧背斜南邊的向斜區,其震源機制解的P軸沿北西西-南東東向,且呈現為走滑型的構造特征(Lei et al., 2017, 2019).由此可見,這兩組位于長寧背斜不同構造部位的地震活動,雖然其震源位置相距僅十幾公里,但中強震震源機制呈現出明顯不同的構造活動特征,表明長寧地區的地震活動具有很強的復雜性(圖1).
震源理論認為,地震是在區域應力場作用下巖體沿某一構造面發生破裂的結果(Scholz, 2002),而震源機制則反映了這一構造活動的具體形態.震源機制解會受到區域應力場的影響,但二者在本質上是不同的,其形態還在很大程度上受到局部發震構造面的影響.那么,長寧地區地震活動的復雜特征是在怎樣的力學成因下形成的呢?僅僅是由于發震斷層構造不同所造成的,還是孕育地震的應力場本身也存在局部變化?這一問題是研究長寧地區地震成因機理的基礎,需要基于區域應力場的可靠結果開展力學分析才能加以解答;而這一應力場結果又需要滿足以下兩方面的條件.首先,區域應力場反演結果需要有足夠的空間精細度,能夠對長寧背斜軸部以及南部向斜區這兩個相距僅十幾公里的不同構造部位之間可能存在的應力場局部變化進行甄別.再者,為了避免應力場反演結果與中強震震源機制解之間的力學對比分析出現自我比較的邏輯問題(例如用中強震震源機制解反演得到的應力場再反過來跟這些中強震震源機制解做力學對比分析),用于反演區域應力場的數據與這些中強震震源機制解應當在數據方面是相對獨立的.目前我們尚未發現滿足上述條件的研究成果見刊,而推進這一研究認識正是本文工作的宗旨所在.

圖1 長寧地區地質構造及中強震震源機制解 黑色圓點表示2014年5月1日至2019年7月10日0.5級以上地震,紅色五角星表示其中5級以上地震,藍色沙灘球表示本文采用gCAP方法反演得到的震源機制解.圖中地震定位結果來自國家地震數據中心的正式觀測報告;地層數據引用自 中華人民共和國地質圖1∶2500000(中國地質調查局, 2004).Fig.1 Geological map and focal mechanism solutions of moderate to strong earthquakes in Changning area Black dots represent the M≥0.5 earthquakes from May 1, 2014 to July 10, 2019, in which red stars mark those M>5.0 earthquakes. Blue beach balls represent the focal mechanism solutions from gCAP inversion in this paper. All the seismic locations in the figure are obtained from the official observation report of the National Seismic Data Center. Stratigraphic data is referred from the geological map of the People′s Republic of China (1∶2500000) (China Geological Survey, 2004).
按照上述分析,本文選定區域應力場的反演方法需要考慮空間精細度以及與中強震震源機制相對獨立這兩方面需求.
利用大量小震P波初動極性求取綜合震源機制解是研究區域應力場的重要途徑(Aki, 1966;許忠淮等, 1987; 萬永革等, 2011a).該方法可以利用大量無法計算出單個震源機制解的小震P波初動極性數據推斷區域應力場的主要特征,因此,其所用數據與中強震震源機制解在數據方面是相對獨立的.同時,小震P波初動極性在數據量和空間分布上優勢顯著,且小震震源尺度小,對局部應力場的變化敏感,使得該方法具有研究應力場精細結構的潛力.鑒于這些原因,本文選用小震綜合震源機制解方法作為開展區域應力場反演的基礎.
然而,如果僅單純地通過小震綜合震源機制解方法進行區域應力場反演,也難以獲取能夠滿足本文區域應力場精細程度的可靠結果,主要體現在以下兩個方面.首先,傳統定位方法給出的小震定位誤差往往有幾公里(張國民等, 2002),這一誤差直接影響著能否正確合理地區分出不同構造部位的小震數據來研究區域應力場的精細變化.再者,地震定位和P波速度結構的誤差還會造成P波初動數據的射線方向誤差,這一射線誤差對于近震數據尤為顯著(Hardebeck and Shearer, 2002);而小震綜合震源機制解方法所用的數據由于震級小,往往正是以近震數據為主(許忠淮等, 1987; 萬永革等, 2011a).雖然有學者認為,在數據十分充足的情況下,數據之間的誤差可能會相互抑制,使得應力場反演結果相對穩定(盛書中等, 2015);但對于精細尺度的應力場反演,需要將反演區域限制在比較小的范圍內,此時P波初動極性數據的充足性難以保證,而數據誤差的影響就難以忽視(胡幸平, 2018).鑒于上述原因,在利用小震綜合震源機制解方法開展區域應力場的精細反演之前,有必要對長寧地區的小震定位以及波速結構進行修正,而雙差層析成像方法(Zhang and Thurber, 2003)能夠很好地契合這一研究需求.該方法是將雙差定位方法(Waldhauser and Ellsworth, 2000)和傳統走時層析成像方法相結合,利用鄰近地震之間相對走時差和地震到臺站之間的絕對走時進行震源位置和波速結構的聯合反演,能夠得到更為準確的地震定位和波速結構(Zhang and Thurber, 2006; 王長在等, 2018).因此,本文將引入雙差層析成像的聯合反演,來改進利用小震綜合震源機制解方法開展的區域應力場精細研究.
本文從國家地震數據中心搜集了2014年5月1日至2019年7月10日發生在27.84°N—28.84°N、114.40°E—115.40°E范圍內0.5級以上地震事件的正式觀測報告,從中提取P波和S波走時,以及P波初動極性.為了保證每個地震事件定位結果的可靠性,我們只選取300 km震中距以內有9條以上走時數據的地震,最終篩選出158個臺站記錄的19886個地震事件的P波走時數據139109條以及S波走時數據137322條(圖2a,b).我們參考Lei等(2017)在長寧地區的研究中所給出的模型設定初始速度模型(圖2c),基于該初始模型利用TauP軟件(Crotwell et al., 1999)計算出的理論時距曲線與篩選出的走時數據具有良好的一致性(圖2b),表明本文的初始模型以及篩選出的走時數據均比較合理.利用這些絕對走時數據,并設定鄰近地震事件對的最大間隔為15 km,我們又得到了P波和S波相對走時差數據分別為808954條和806272條.
基于上述絕對走時以及相對走時差數據,本文利用Zhang和Thurber(2003)開發的tomoDD程序對長寧地區的地殼波速結構和震源位置進行了聯合反演.反演采用帶阻尼的LQSR算法,阻尼因子和光滑權重因子通過繪制數據方差和模型方差的均衡曲線來確定其最優值(Aster et al., 2012; 王長在等, 2018),分別為500和40.此外,通過棋盤格檢測板測試(Humphreys and Clayton, 1988),我們發現長寧背斜區波速結構的橫向反演分辨率能夠達到0.1°×0.1°,因此將這一核心區域的反演分辨率設定為0.1°×0.1°(圖3).
通過反演,本文得到了15762個地震的重定位結果(圖4),其橫向及深度的平均相對誤差分別為422 m和534 m.相比于初始定位結果(圖1),重定位結果在平面上更為集中,長寧MS6.0地震及其5級以上中強余震更清晰地呈現出沿西北方向線性分布的特征.在震源深度上,重定位后的地震主要集中在2~6 km,與Meng等(2019)對2015年2月至2017年12月發生在28°N—28.25°N、104.6°E—105°E范圍內的地震所做的定位結果在總體上是一致的.如果進一步區分震級,本文結果還顯示出M<4的小震震源深度仍集中在2~6 km,但M≥4的中強震則更多地分布在6~9 km.對照何登發等(2019)基于反射地震剖面給出的地層結構,可以看出,前者主要分布在基底之上的震旦系及以上地層,而后者主要分布在震旦系以下的基底中.

圖2 (a)地震射線分布圖; (b)震相走時曲線; (c)初始速度模型Fig.2 (a) Distribution map of seismic rays; (b) Travel-time curve; (c) Initial velocity model

圖3 長寧背斜區棋盤格檢測板測試結果(預設的速度擾動幅度為±5%)Fig.3 Results of checkerboard test with ±5% velocity perturbation around Changning anticline
與此同時,反演也給出了長寧地區的地殼波速結構.基于檢測板測試的恢復分辨率大于0.6的反演結果,本文通過插值得到了長寧背斜區0.01°×0.01°的地殼波速結構(圖5).由這一結果我們注意到,長寧背斜區在7 km以上主要表現為明顯的S波高速異常.對此,本文整理前人對長寧地區地層構造所做的研究發現,長寧背斜區的地層厚度可達9 km(郭正吾等, 1996; 何登發等, 2019);在漫長的地質時期中,由于受到來自西南方向的擠壓應力,長寧背斜區的地層發生了褶皺、抬升,后又經歷地表的風化剝蝕,形成了長寧背斜現今的構造面貌(何登發等, 2019).據此本文推測,正是由于地層的前期抬升和后期剝蝕,造成了長寧背斜區7 km以上的地層與周邊區域更深部的地層相對應;由于S波速度一般會隨著深度增加,那么這種橫向上的地層差異就導致了在同一深度上長寧背斜區的S波速度相比于周邊區域呈現為相對的高速異常.雖說由于鹽氣開采長寧地區存在大量的注水現象(Sun et al., 2017; Lei et al., 2019),可能會造成S波低速異常(Roland et al., 2012; Guo et al., 2018),但從本文結果來看,這些注水作用對于當地S波速度結構的影響并不顯著,并沒有顛覆不同地層S波速度的高低差異.

圖4 長寧地區地震重定位平面(a)及深度分布(b)Fig.4 Horizontal (a) and vertical (b) distribution of the relocated earthquakes in Changning area

圖5 不同深度反演結果與初始模型的相對差異 黑色圓點表示臨近深度內的地震重定位結果,紅色五角星表示其中5級以上地震.Fig.5 Relative differences between inversion result and initial model at different depths Black dots indicate the relocated epicenters close to each depth, in which red stars indicate M>5.0 earthquakes.
此外,我們也嘗試了其他的初始速度模型(趙珠等, 1997)進行反演,得到了非常相似的地震定位和波速結構的反演結果,表明本文的雙差層析成像反演結果對于不同的初始速度模型來說是比較穩定的.
在上述反演結果的基礎上,本文僅選取長寧MS6.0地震發生前的M<4.0小震初動極性數據進行應力場反演,這一方面是為了保證應力場反演所用的數據與中強地震震源機制解的獨立性,另一方面則是為了排除長寧MS6.0地震可能造成的局部應力調整.依據雙差層析成像得到的P波速度結構及地震重定位結果,本文采用偽彎曲法(Um and Thurber, 1987)計算了這些極性數據的射線方向,獲得了6936條修正后的P波初動極性數據.
本文首先采用網格掃描的方法對區域應力場的橫向不均勻性進行了掃描.我們將長寧背斜區劃分為0.1°×0.1°的網格;對于每一網格點,出于對應力場反演的區域平滑度和局部保真度的雙重考慮,選取橫向距離(DL)小于20 km的極性數據,并對其設定三次衰減函數的距離權重因子(WD):
WD=(1-(DL/DL0)3)3,
(1)
式中DL0為橫向距離的閾值,等于20 km.此外,考慮到清晰和不清晰的初動極性通常具有不同的拾取誤差(Hardebeck and Shearer, 2002),我們設定其清晰度權重因子分別為1和0.5.本文采用格點嘗試方法(許忠淮等, 1983; 俞春泉等, 2009)來求取小震綜合震源機制解,搜索步長為2°×2°×2°,選取加權矛盾比最小的200個解作為可選解.通過對這些可選解的三個主應力軸構成的張量進行張量平均,然后求解平均張量的主向量,本文最終確定出網格點上每個主應力軸的平均方向;此外,通過計算網格點上每個主應力軸與可選解的對應主應力軸之間的空間夾角的均方根,本文得出網格點上每個主應力軸的方向離散度.為了保證結果可靠性,本文只對極性數據超過500條的應力反演結果進行展示和分析(圖6).
從圖6可以看出,除了少數邊緣網格點之外,大部分結果的主應力軸的方向離散度都在15°以內,表明這些反演結果是相當穩定的.區域應力場的網格掃描結果顯示出,長寧地區地殼應力場的最大主應力軸(σ1)在整個研究區內基本都呈現為低傾角的近水平狀態,其方位大體沿著近東西向,但在長寧背斜軸部往北以及南部向斜區這兩個區域之間存在不太明顯的分叉現象——前者偏向北東東-南西西向,而后者偏向北西西-南東東向.相比之下,這兩個局部區域應力場在構造應力類型上的差異更為顯著,其中長寧背斜軸部主要呈現為中等主應力軸(σ2)低傾角、最小主應力軸(σ3)高傾角的逆沖型應力結構,而南邊的向斜區則主要呈現為σ2高傾角、σ3低傾角的走滑型應力結構.
上述局部應力差異大致以28.3°N的緯度線為界,為了更清晰地對比這兩個區域的應力場差異,本文將修正后的P波初動極性數據分成28.3°N以南和以北兩組(簡稱為南組和北組),對于每個數據組內不考慮距離權重,同樣采用上述格點嘗試法和參數設定求解小震綜合震源機制解(圖7).分組結果顯示出,南北兩組的最大主應力軸(σ1)的平均傾角都較低(南組10.1°,北組12.9°),接近水平狀態,其平均方位雖然從北向南顯示有17°左右的順時針偏轉(北組79.4°,南組276.6°),但從二者可能解的σ1方位分布來看,兩組σ1方位差別并不顯著,可以說都大致沿著近東西向.與σ1軸相比,σ2和σ3軸呈現出的分組差異明顯超出了其可選解的分布范圍,其中北組表現為σ2低傾角、σ3高傾角的逆沖型應力結構,而南組則表現為σ2高傾角、σ3低傾角的走滑型應力結構,而且這種構造應力類型的差異是相當顯著的.
本文反演得到的南部向斜區內的最大主應力方位與Lei等(2017, 2019)利用該地區震源機制解的反演結果以及王玉滿等(2016)整理長寧氣田區地應力測量資料得到的最大水平主應力方位比較一致,而長寧背斜軸部的最大主應力方位則與Zheng等(2018)和石玉濤等(2013)給出的該地區北東-南西向快波方向比較接近.

圖6 長寧地區網格點上主應力方向及其離散度 黑線代表主應力方向,其長度與主應力軸的傾角成反比,色階表示每個主應力軸的方向離散度. 黑色圓點表示提供P波初動極性數據的M<4.0小震重定位結果.Fig.6 Orientations of principal stress axes and their direction dispersions determined by grid-scanning in Changning area Black bars indicate the orientations of principal stress axes, and their lengths are inversely proportional with the plunges of the principal stress axes. The direction dispersion of each principal stress axis is color coded. Black dots indicate the relocated M<4.0 earthquakes providing initial motion data of P wave.

圖7 長寧地區應力場分組反演結果 左圖中圓點表示提供P波初動極性數據的M<4.0小震重定位結果;右圖中-、o、+分別表示可選解的σ1、σ2、σ3軸,平均解的 σ1、σ2、σ3軸的方位角和傾角標注在對應應力軸符號之下.北組和南組的地震和應力場反演結果由紅色和綠色加以區分.Fig.7 Stress inversion results determined by discriminating seismicity groups in Changning area Dots in the left map indicate the relocated M<4.0 earthquakes providing initial motion data of P wave. -, o, + in the right polar graph indicate optional σ1,σ2, and σ3 axes respectively; the average σ1,σ2, and σ3 axes are marked namely by the symbols of σ1,σ2, and σ3, and their azimuth and plunge are written below. Red/green color discriminate the relocated earthquakes and stress inversion results of north/south seismicity group.
為了研究區域應力場與中強地震的力學一致性,同時也為了進一步驗證局部應力差異,本文將上述兩個局部區域的應力結構與中強震震源機制解進行了力學對比分析.對此,本文采用的是計算歸一化剪滑分量(Angelier, 2002)的方法.歸一化剪滑分量(ω)同時考慮兩個參數,一是相對剪應力大小,即區域應力場(應力張量)在震源機制解節面上的計算剪應力除以最大剪應力(胡幸平, 2018; 萬永革, 2020),二是剪滑角,即實際震源機制解節面上的滑動方向和計算得到的剪應力方向之間的角度(Gephart and Forsyth, 1984; 許忠淮, 1985; Michael, 1987).歸一化剪滑分量(ω)的變化范圍為-1~1,值越小,表明很難用區域應力場解釋實際震源機制解所指示的巖體滑動,若為負值,則表明區域應力場與實際震源機制解在力學上是相抵觸的.本文之所以使用歸一化剪滑分量來衡量區域應力場與震源機制解的一致性,而不是單純使用剪滑角,是因為后者僅考慮了如果巖石在區域應力場下發生滑動,那么這種滑動與實際震源機制解之間的差異,而前者還同時從剪應力量值的角度考慮了巖石在區域應力場下發生滑動的能力或者說概率.如果區域應力場作用在震源機制解節面上的相對剪應力很小,說明巖石滑動很難被驅動,那么即便是剪滑角很小,也難以認為區域應力場與震源機制解在力學上是一致的.另外,歸一化剪滑分量的計算對于震源機制解兩個節面中的任意一個,其結果都是相等的,這也避免了從震源機制解兩個節面中判斷發震斷層面的難題.
本文用于對比分析的中強震震源機制解是我們采用國際上常用的gCAP方法(Zhu and Ben-Zion, 2013)反演得到的,所用的波形數據來自中國地震局地球物理研究所國家測震臺網數據備份中心.在震源機制反演中,Pnl和S波截取波形窗長分別為35 s和70 s,相應濾波范圍為0.02~0.2 Hz和0.02~0.1 Hz;走向、傾角和滑動角的搜索步長為10°,深度為1 km.格林函數采用頻率-波數法(FK)來計算(Zhu and Rivera, 2002),其中采樣間隔設為0.1 s,采樣點為1024個.震源機制解的反演結果見圖1及表1.

表1 長寧地區中強震震源機制解反演結果Table 1 Focal mechanism solutions of moderate to strong earthquakes in Changning area
計算歸一化剪滑分量需要用到相對應力大小R;按照Gephart和Forsyth(1984)以及萬永革等(2011b)的定義,其表達式為
R=(σ2-σ3)/(σ1-σ3),
(2)
R的變化范圍為0~1,但利用小震綜合解的反演并不能確定出其具體參數(許忠淮等, 1983; 萬永革等, 2011b).對此,本文首先將其設定為0.5這一中間值來進行計算分析(圖8).結果顯示出,對于對應區域內的應力場和中強震震源機制解(例如同樣是長寧背斜軸部的應力場反演結果和中強震震源機制解),計算得到的歸一化剪滑分量都比較大(13個震源機制解中有10個的ω值大于0.6,ω最小值為0.29),這就表明二者之間在力學上具有良好的一致性,同時這也意味著驅動小震和中強震的應力場是幾乎同一的;與之相比,對于非對應區域的應力場和中強震震源機制解(例如長寧背斜軸部的應力場反演結果和南部向斜區內的中強震震源機制解),計算得到的歸一化剪滑分量明顯偏小(13個震源機制解中有12個的ω值小于0.3,甚至有6個ω值為負值),這就表明二者之間在力學上的吻合度較低,甚至是相抵觸的.這一強烈的對比結果表明,長寧背斜區的地殼應力場確實存在比較明顯的局部差異,而這種應力場的局部改變是造成背斜軸部和南部向斜區內的中強震活動機制存在差異的必要力學基礎.

圖8 長寧地區中強震震源機制解與 區域應力場的力學一致性Fig.8 Results of mechanical fitting between the stress fields and focal mechanism solutions of moderate to strong earthquakes in Changning area
從圖8中也能發現第11號震源機制解在非對應計算中反而會得到更大的歸一化剪滑分量.本文認為此個別力學擬合結果的反置可能由于震源機制解反演誤差或者其反映的是應力場調整等原因導致的,而這種極個別現象并不能顛覆應力場與中強震震源機制的對比分析認識.
對于長寧地區的相對應力大小R,Lei等(2019)以及劉敬光等(2019)通過震源機制解反演給出了相對較低的取值:前者結果為0.15;后者在R=(σ2-σ1)/(σ3-σ1)的定義下的結果為0.7,如果按照本文R的定義,其結果為0.3.因此,本文也將R設為0.15、0.3等低值,以及試算了0.8這一高值,此外還采用其他學者的震源機制結果(易桂喜等, 2019)進行了多組計算分析.這些設定不同相對應力大小以及采用不同來源的震源機制數據的計算結果都表現出統一的對比趨勢,即利用小震綜合震源機制解反演得到的應力場和中強震震源機制解在對應區域內一致性較好,但在非對應區域之間一致性明顯降低.因此,可以認為本文的對比分析是相當穩定的,而且這一力學分析能夠說明長寧背斜軸部和南部向斜區之間確實存在局部應力場的改變.
通過長寧背斜區應力場的精細反演及其與中強震震源機制解的力學對比分析,本文發現,長寧地區地殼應力場存在顯著的局部變化,長寧背斜軸部表現為最大主應力軸北東東-南西西向的逆沖型應力結構,而南部向斜區則表現為最大主應力軸北西西-南東東向的走滑型應力結構.那么在相距僅十幾公里的兩個構造部位,應力場為何會表現出這樣的差別?
細致觀察上述局部應力結構(圖7)可以看出,雖然二者在最大主應力方位上也有一定的差異,但更顯著的差異在于前者的中間主應力軸近水平、最小主應力軸近直立,而后者與之相反.這就表明相比于南邊的向斜區,長寧背斜軸部在橫向上受到的擠壓作用更強.本文認為這種橫向擠壓作用的差異可能是由于巖體重力在橫向上產生的圍壓不同而帶來的.根據何登發等(2019)對長寧背斜區地層系統的研究,長寧背斜軸部10 km以上的地殼主要由寒武系、震旦系以及更老的結晶基底構成,巖體中包含相當比例的震旦系白云巖和基底變質雜巖;而南部向斜區10 km以上的地殼巖層中有3 km以上厚度的奧陶系至侏羅系沉積巖,巖層組成中頁巖所占比例高于前者,而也正是如此,南部向斜區成為長寧氣田的集中開采區(王玉滿等, 2016).在巖石的力學性質上,一般頁巖的泊松比明顯低于白云巖和變質雜巖,因此有理由認為,長寧背斜軸部巖體的重力由于泊松效應會在橫向上產生更大的圍壓,而這一橫向圍壓疊加在周邊構造活動產生的橫向推擠作用后,使得水平方向上的最大和最小主應力(SH、Sh)增加從而成為最大及中等主應力(σ1、σ2),而垂向應力(SV,可能主要由上覆巖體重力產生)成為最小主應力(σ3);相比之下,在南部向斜區,巖體重力由于泊松效應產生的橫向圍壓較前者明顯更小,垂向應力(SV)仍為三個主應力中的中等主應力(σ2).于是這種由于巖石泊松比差異帶來的橫向圍壓差異就導致了長寧背斜軸部和南部向斜區這兩個局部區域之間存在逆沖型(SH>Sh>SV)和走滑型(SH>SV>Sh)的應力類型轉變.
為了檢驗上述成因的可能性,本文依據雙差層析成像反演的波速結構計算波速比(VP/VS),然后依據波速比與泊松比(μ)的關系,即
μ=(0.5×(VP/VS)2-1)/((VP/VS)2-1),
(3)
計算了長寧背斜區的地殼泊松比分布(圖9).
從圖9可以看出,長寧背斜區的巖石泊松比具有明顯的橫向差異;這種差異在7~10 km的深度上表現得尤為顯著,長寧背斜軸部的泊松比明顯高于南部向斜區,差值甚至可以達到0.15左右,而這種巖石力學性質的強烈差異很可能就是由地層的橫向差異所造成的.另一方面,王玉滿等(2016)通過整理長寧氣田的地應力數據指出,氣田區2300~3200 m深度上的橫向差應力值(SH-Sh)為21.4~22.3 MPa,垂向應力(SV)為56~66 MPa.依據上述這些數據,本文對長寧背斜區的地應力量值以及橫向圍壓差異的影響進行估算.
吳珍漢和白加啟(1997)根據多種應力資料的統計分析指出,上地殼的最大差應力在25 MPa左右;由此看來,長寧氣田區深度3 km左右的橫向差應力值(22 MPa)已經處于超高壓狀態.因此,有理由認為該地區上地殼內的橫向差應力不會再隨深度明顯增加,于是我們將深度7 km左右的橫向差應力按25 MPa估計.選擇長寧氣田所在的南部向斜區的應力結構作為區域應力場的參考基準,根據反演結果,該區應力結構呈現為走滑型,那么這種橫向差應力(SH-Sh)與(σ1-σ3)接近.垂向應力(SV)在深度3 km左右約為60 MPa,與上覆巖體的重力相當,據此,我們將深度7 km左右的垂向應力(SV)估算為140 MPa.背斜軸部和南部向斜區的泊松比差值按0.1估算,那么前者由于泊松效應在橫向上產生的圍壓比后者高出14 MPa左右.如果相對應力大小R按0.5估算,則σ2-σ3=(σ1-σ3)×R≈12.5 MPa,可以看出這一由于泊松比不同導致的~14 MPa的橫向圍壓差別具備了造成σ2和σ3在Sh和SV之間轉換的量值需求;而如果相對應力大小按照0.15(Lei et al., 2019)或者0.3(劉敬光等, 2019)來估算,則σ2-σ3等于3.75 MPa或者7.5 MPa,那么橫向圍壓差別更是明顯超出該值,就更有可能導致σ2和σ3在Sh和SV之間的轉換.

圖9 長寧背斜區波速比及泊松比分布 黑色圓點表示臨近深度內的地震重定位結果,紅色五角星表示其中5級以上地震.Fig.9 Poisson′s ratio and velocity ratio at different depths around Changning anticline Black dots indicate the relocated epicenters close to each depth, in which red stars indicate M>5.0 earthquakes.
通過上述應力量值的估算,可以看出,背斜軸部和向斜的巖體在泊松比這一力學性質上的差異,確實能夠造成這兩個局部應力場具有不同的應力類型.同時,這種巖石泊松比差異造成的圍壓差在長寧背斜軸部到向斜之間的橫向作用,很難說能夠與區域應力場在水平面上的任一主應力(SH或Sh)方向一致,于是也會導致局部最大主應力方向發生一定程度的改變(Zoback,1992).前人的應力場數值模擬研究也證明了巖石介質差異具備造成應力類型和應力方向發生改變的這種可能性(劉平江等, 2007).
另外,結合圖6和圖9還可以發現,南部向斜區內部也依然存在更小尺度的局部應力結構和巖石泊松比的對應差異.例如,興文MS5.7地震和珙縣MS5.1地震附近呈現為走滑型應力結構和低泊松比,而其西南邊(104.7°E、28.1°N附近)呈現為逆沖型應力結構和高泊松比.由此可見,在長寧背斜軸部與南部向斜區之間、以及南部向斜區內的更小尺度上,局部應力結構和巖石泊松比之間都表現出很強的對應關系,這就進一步支持巖石泊松比的橫向差異很可能就是局部應力結構差異的主要成因.
除了巖石泊松比的橫向差異外,其他更為復雜的因素也有可能在一定程度上引起或者說加劇了長寧背斜區的應力場局部變化.例如,長寧背斜區由于鹽氣開采不斷發生的高壓流體注入.鑒于背斜軸部采鹽和南部向斜區采氣的注水工藝可能不同、流體滲透可能受裂隙發育方向的影響、以及注水過程的持續動態性,該地區的流體注入過程造成的地應力變化很可能不是靜水壓力那樣的各向同性,就有可能導致這兩個區域的應力方向和結構形態產生差異.再如,現今構造作用以及巖體中殘留的古應力場釋放在長寧背斜這一古褶皺構造的不同部位是否也存在局部變化(常祖峰等, 2020),從而導致了現今應力結構的不同,而這種因素的影響也不能完全排除.
上述這些都是造成區域應力場局部改變的可能因素;同時,另一方面,無論是用小震極性還是震源機制解反演得到的都是地震釋放的應力,其主應力軸會受局部斷層構造影響,可能與實際應力場的主應力軸有所差別(Tian et al., 2013).雖說通過大量不同走向和傾角的斷層上的地震數據能夠有效地推斷出區域應力場的方向特征(許忠淮等, 1983),且本文所用的小震數據豐富,可以認為上述條件大致成立,但仍難以保證本文反演結果是對實際應力場的“完美”恢復.這種地震釋放應力與實際應力場的可能差異也是目前所有基于地震學的應力場反演方法無法完全消除的問題,而這也有可能是造成應力場反演結果出現南北差異的部分原因.
通過上述分析可以得出,巖石泊松比的橫向差異具備造成長寧背斜區地殼應力場的局部變化的可能性,因此,這種局部應力變化是可以理解的.但與此同時,這種局部應力變化也有可能是多方面的復雜因素共同作用的產物,甚至還可能包含了少量由于地震學反演造成的“假象”,因此,需要我們結合多種數據和手段、開展更為全面和深入的研究工作才能徹底解答.
基于上述雙差層析成像聯合反演結果、震源機制解、以及區域應力場結果,本文也對長寧地震的發震構造進行了探討.
從重定位和震源機制解的結果可以發現,長寧地震序列的中強震沿著北西方向線性分布,其震源機制解也幾乎都具有北西向的節面,因此,本文選取了2條穿過長寧MS6.0地震的剖面來討論這一地震的發震斷層構造.其中一條剖面沿著長寧地震序列線性分布的北西方向(圖4中的AA′,走向NE130°),另一條與之正交(圖4中的BB′,走向NE40°);開展剖面分析的底圖數據包括2個,一是地震相對密度(圖10),另一則是波速結構(圖11).
本文地震相對密度的計算方法為:將整個研究區劃分為0.01°(經度)×0.01°(緯度)×1 km(深度)的三維網格,對于每個網格點,統計其周邊2 km范圍內不同震級區間的地震數量,然后通過除以所有網格點周圍對應震級區間地震數量的最大值進行歸一化,得出各個網格點上的不同震級區間的地震相對密度. 從圖10的AA′剖面可以看出,震源深度整體上向北西方向逐漸加深,且不同震級地震的總體震源深度存在分離現象,中強震主要分布在小震底部.對照何登發等(2019)給出的地層結果,可以看出4級以上的中強震主要分布在基底中.由于這些中強震的破裂長度可達數公里(Wells and Coppersmith,1994),且它們的發震節面和空間位置具有一定的連貫性,因此,本文認為在長寧背斜軸部可能存在北西走向的基底斷層,控制著4級以上中強震的活動.這種基底斷層可能是一條較大規模(長度超過15 km)的斷層,又或是多條中等規模(長度數公里)的斷層組合而成;但無論哪種方式,其總體上都基本貫通(或者說組合貫通)了長寧背斜的軸部.依據AA′剖面上的地震相對密度和中強震震源機制解(圖10),并結合該剖面的S波速度結構(圖11),本文勾畫了這種基底斷層沿著北西向AA′剖面的構造特征,認為其主要分布在6~9 km深度內,且從長寧MS6.0地震處往西北方向略微加深.相比之下,小震幾乎都分布在中強震震源以上的地殼淺部,其優勢分布深度與中強震存在一定程度的分離,且小震破裂尺度小,因此本文認為這些淺部小震可能并不對應著規模較大的基底斷層的活動,而是上覆地層中的裂隙或者微小斷層活動而產生的.

圖10 垂直剖面上的地震相對密度及震源機制投影 AA′和BB′剖面位置見圖4;紅色虛線表示推斷的斷層構造.Fig.10 Relative density of seismicity and projection of focal mechanism solutions on vertical profiles The locations of profiles AA′ and BB′ are sketched in Fig.4. Dashed red lines indicate the inferred seismogenic fault.

圖11 垂直剖面上的波速結構及震源機制投影 AA′和BB′剖面位置見圖4;黑點表示剖面兩側2 km以內的地震,紅點表示4級以上地震;紅色虛線表示推斷的斷層構造.Fig.11 Seismic velocity and projection of focal mechanism solutions on vertical profiles The locations of profiles AA′ and BB′ are sketched in Fig.4. Black dots indicate the relocated earthquakes within 2 km on both sides of each vertical profile, while larger red dots indicates those M≥4.0 events. Dashed red lines indicate the inferred seismogenic fault.
另外,BB′剖面的波速結構似乎呈現出向西南緩傾的低波速異常(圖11),有可能指示著斷層附近的巖體破碎.因此,本文猜想長寧背斜軸部的基底斷層可能是傾向西南,傾角較緩,在40°左右,而這種傾角的推測也與本文以及易桂喜等(2019)的震源機制解節面的統計特征基本相符.
基于上述斷層構造推測(走向130°,傾角40°),利用長寧背斜軸部的應力場反演結果,分別將相對應力大小R取為0.15、0.3、0.5以及0.8,可以計算出該斷層構造上的理論滑動角和相對剪應力大小(胡幸平, 2018; 萬永革, 2020).這些計算結果(表2)表明,對于上述R的任一不同取值,計算出的斷層面上的相對剪應力大小都超過0.9,表明本文推測的斷層構造在區域應力場下具有很強的錯動能力.另外,隨著R取值的增加,計算出的斷層錯動的逆沖分量逐漸增加;對于R≤0.5的三個取值,其斷層錯動性質都是左旋走滑與逆沖比例相當,與實際長寧地震序列中強震震源機制解的總體特征相當吻合;而即使是R取為0.8,計算出的斷層錯動性質與實際情況相比也沒有本質差別.因此,本文認為長寧地震的發生正是由于背斜軸部存在的基底斷層在局部區域應力場作用下發生錯動的結果.

表2 采用不同相對應力大小計算出的推斷斷層構造面上的 理論滑動角和相對剪應力大小Table 2 Theoretical rakes and relative shear stresses calculated by using different relative stress ratios on the inferred fault structure
本文通過雙差層析成像聯合反演獲取了長寧地區的小震定位和波速結構,在此基礎上利用小震綜合震源機制解方法反演了長寧地區地殼應力場的精細結構.結果顯示出,長寧地區地殼應力場的最大主應力軸在整個研究區內基本都呈現為近水平狀態,其方位雖然由北向南發生了一定角度的順時針旋轉,但也基本保持在近東西向;構造應力類型在長寧背斜軸部和南部向斜區之間存在顯著差異,前者表現為逆沖型應力結構,而后者表現為走滑型應力結構.
通過上述應力場結果與gCAP反演得到的中強震震源機制解做力學對比分析,本文發現長寧背斜軸部與南部向斜區的局部應力結構與其各自對應區域內的中強震震源機制解的吻合度較高,但與非對應區域內的中強震震源機制解的吻合度較低,甚至在力學上是相抵觸的.據此本文認為,區域應力場的局部改變是長寧地區中強震震源機制在十幾公里范圍內就發生顯著變化的必要力學基礎.而通過巖石力學估算分析,本文認為長寧地區的局部應力變化很可能主要是由巖石泊松比的橫向差異造成的,但也有可能是多方面因素的共同產物.
最后,依據小震精定位和波速結構的剖面投影以及基于區域應力場的估算分析,本文推斷了長寧地震的發震構造,認為長寧背斜軸部可能存在基底斷層,控制著包括長寧MS6.0地震在內的中強地震的活動;而長寧地區地殼淺部的小震活動與中強震在深度上存在比較明顯的分離,可能是由于上覆地層中的裂隙或者微小斷層活動而產生的,而與軸部基底的較大規模的斷層的關系較弱.
本文對長寧地區的地殼應力場及地震活動的研究分析側重在空間差異上,而在不同的時間段,例如,在鹽氣開采的不同時期、以及長寧MS6.0地震等重要地震事件前后,是否也存在變化,由于資料所限并沒有展開研究.然而,這種時間上的變化對于該地區地震成因機理研究和活動趨勢預判都是相當重要的,值得繼續深入研究.
致謝本文研究中所用的地震觀測報告來自國家地震數據中心,波形數據來自中國地震局地球物理研究所國家測震臺網數據備份中心,繪圖通過GMT軟件包(Wessel and Luis, 2017)完成;應急管理部國家自然災害防治研究院的徐錫偉研究員和雷建設研究員對本文工作提供指導意見;匿名專家對本文提出修改意見.在此感謝上述機構和個人對本文工作的支持.