日本至中國東北區域構造應力場分布與Benioff帶的形態關系研究

馮 兵 李 楊 王文濤 楊應召 汶宇龍 劉 蜀

1 中國地震局第二監測中心，西安市西影路316號，710054

自板塊構造學說提出以來，環太平洋地震帶一直是地震學家和地質學家研究的重點地區。對環太平洋俯沖帶的研究是認識板塊構造運動的有力證據，而發生在俯沖帶的地震則是認識板塊構造運動及板塊間應力場狀態的最有效手段。地震的本質是板塊運動跨地質斷層的錯動，從純物理的角度來看，地震是一種非平衡系統的耗散，地震的發生是一個能量釋放的過程。

日本海溝是太平洋板塊向歐亞板塊俯沖，同時受菲律賓板塊和歐美板塊共同作用的結果(圖1)，是全球地震活動非常頻繁的區域之一。海溝俯沖帶的發現是板塊構造運動的有力證據。Wadati[1]研究發現，發生在大洋中的深源地震會形成按一定規律排列的地震區，在此研究基礎上發現了著名的Wadati-Benioff帶。寧杰遠等[2]利用1965～1982年的地震數據資料分析日本至中國東北地區地震深度的分布，認為日本本州至中國東北Benioff帶具有連續性，中國東北區域深震為日本海下Benioff帶的一部分。

在日本3·11 MW9.0大地震發生以后，國內外許多學者對該區域的板塊應力狀態及同震效應進行了詳細研究，并對日本海溝俯沖帶區域的應力場變化進行分析，以期解釋地震發生的相關物理過程。Yokota等[3]利用多種數據資料反演3·11大地震的破裂過程，結果顯示，地震破裂長度為480 km，寬度為180 km，同震最大位移量達到35 m。陳為濤等[4]研究認為，該地震的同震位移在瞬間產生的形變量相當于中國東北地區12.7 a的長期應變積累。在俯沖帶區域，主壓應力軸在大震之前處于板塊匯聚方向，大震之后地殼上盤及部分上隆區域出現張應力變化[5]。楊佳佳等[6]通過震源機制解反演日本3·11大地震震源區應力場變化發現，震源區主壓應力軸和最大水平主壓應力都發生了改變。

上述研究大多通過劃分地震的區域和時間進行構造應力場分析，但關于Benioff帶的空間三維構造應力場分析及與中國東北地區構造應力場的相關研究較少，因此難以從深度及空間方面分析該地區應力場變化的整體分布。為了更加科學和詳細地認識日本海溝俯沖帶的應力場變化及中國板塊應力轉換的過程，本文將基于震源機制解的數據資料進行空間應力場3D反演。通過分區域、分空間進行應力場反演，得到日本至中國東北地區的構造應力場及其空間變化，并以此為基礎，探討海溝俯沖帶的具體分布情況，同時分析日本海溝俯沖帶與中國東北Benioff帶地下構造的相互作用，以期從更大范圍認識板塊運動的變化情況，為地震預測預報提供參考。

黑色實線為板塊邊界，黑色虛線表示正在逐漸成形的塊體邊界，藍色、紅色圓圈及綠色五角星為震中分布，矩形框為研究區圖1 日本海溝至中國東北區域地形構造及研究區震中分布Fig.1 Topography structure and epicenters distribution of the study area from Japan trench to northeast China

1 數據資料與研究方法

1.1 數據資料選取

本研究收集了1997-01-01～2020-05-07MW≥3.0共2 291個地震的震源機制解數據資料，所選數據緯度為35～50°N，經度為122～140°E，深度在0～700 km之間。所有數據資料均來源于NIED(National Research Institute for Earth Science and Disaster Prevention) F-net寬頻帶地震臺網中心矩張量解目錄。

1.2 研究方法

應用震源機制解推斷某一地區的應力場方向，是認識地殼運動和地震震源物理過程的有效方法之一。SATSI程序采用阻尼最小二乘法，可同時得到多個網格點的應力張量并進行平滑，具有良好的計算性能[7]。其基本原理為：最小化每一個斷層面上的滑動矢量與剪應力的差別，建立斷層面震源機制解參數與應力張量的線性關系：

Gm=d

式中，矩陣G可由每個震源機制解的斷層面滑動矢量得到，m為應力張量的單位模型矩陣向量，d為數據矩陣矢量。通過最小二乘法進行應力反演，其中應力張量的模型矩陣向量可確定最大主應力、中間主應力、最小主應力的方向及表示主應力相對大小的R值。本研究采用基于SATSI法發展而來的MSATSI構造應力反演方法[8]，可在MATLAB環境下進行完整的應力反演。由于本次研究基于3D反演進行區域應力場反演，每個網格點與周邊6個網格點相鄰，采用阻尼擬合參數會對每個網格點造成較大影響。應力場反演可得到反映應力相對大小的R值的量θ(θ=1-R=(σ2-σ3)/(σ1-σ3))和3個主應力的方位角與傾伏角，同時還可反演每個網格的極射赤平投影P、T軸分布圖。R取不同值，相應的應力主軸表現結果會存在差異。R值可用來評估主應力軸方向的不確定性，取值范圍為0～1。當R取值為0時，即主張應力軸方向確定，主壓應力軸(σ1)取向與中間應力軸(σ2)取向趨于一致，并在垂直于主張應力軸(σ3)的平面內呈現擠壓狀態；當R取值為0.5時，表明主張應力軸與主壓應力軸確定，相應的應力軸分辨狀態趨于明顯；當R取值為1時，表明主壓應力軸方向確定，主張應力軸取向與中間應力軸取向趨于一致，并在垂直于主壓應力軸的平面內旋轉，同時呈現張性特征；R值在0～0.5時，中間應力軸的壓性特征越來越明顯；R值在0.5～1時，中間應力軸的張性特征越來越明顯[6,9-10]。

2 數據整理與分析

2.1 數據資料分區說明

為從空間上給出日本至中國東北區域的構造應力場，以經度、緯度和深度為劃分標準，對收集到的震源機制數據進行分區，將分區后的震源機制解數據置于3×4×5的矩陣網格中。將相應的網格矩陣單元作為固定區域的應力場反演區塊，以序號為i、j、k的網格點對反演的區域進行命名分區。對于反演坐標軸中的X軸與Y軸，以5°×5°的經緯度進行區域劃分。緯度在35°～40°劃分為矩陣[i,0,k]區(i、j、k表示對應矩陣網格的未知數)，緯度在40°～45°劃分為矩陣[i,1,k]區，緯度在45°～50°劃分為矩陣[i,2,k]區。經度小于125°劃分為矩陣[0,j,k]區，經度在125°～130°劃分為矩陣[1,j,k]區，經度在130°～135°劃分為矩陣[2,j,k]區，經度在135°～140°劃分為矩陣[3,j,k]區。深度在0～30 km劃分為矩陣[i,j,5]區，深度在30～60 km劃分為矩陣[i,j,4]區，深度在60～100 km劃分為矩陣[i,j,3]區，深度在100～300 km劃分為矩陣[i,j,2]區，深度在300～500 km劃分為矩陣[i,j,1]區，深度在500～700 km劃分為矩陣[i,j,0]區，具體分區見表1。雖然部分矩陣區塊可能數據量小或缺少震源機制解數據資料，但不會對該區域整體應力場反演的過程及結果造成影響。根據矩陣分區規則，對各個存在震源機制解的網格區域進行統計及劃分：矩陣[0,0,5]區分配到1個震源機制解；矩陣[1,0,5]區分配到28個震源機制解；矩陣[2,0,5]區分配到240個震源機制解；矩陣[3,0,5]區分配到1 083個震源機制解；矩陣[3,1,5]區分配到137個震源機制解；矩陣[3,2,5]區分配到2個震源機制解；矩陣[1,0,4]區分配到2個震源機制解；矩陣[3,0,4]區分配到276個震源機制解；矩陣[2,1,4]區分配到1個震源機制解；矩陣[3,1,4]區分配到1個震源機制解；矩陣[3,0,3]區分配到111個震源機制解；矩陣[3,0,2]區分配到169個震源機制解；矩陣[3,1,2]區分配到75個震源機制解；矩陣[2,0,1]區分配到30個震源機制解；矩陣[3,0,1]區分配到75個震源機制解；矩陣[2,1,1]區分配到7個震源機制解；矩陣[3,1,1]區分配到23個震源機制解；矩陣[3,2,1]區分配到3個震源機制解；矩陣[2,0,0]區分配到4個震源機制解；矩陣[2,1,0]區分配到20個震源機制解；矩陣[2,2,0]區分配到1個震源機制解；矩陣[3,2,0]區分配到2個震源機制解。為保證數據反演結果的可靠性，設置約束反演應力場的相關條件，見表2。

表1 矩陣網格劃分標準

表2 應力反演的相關約束條件

在構造應力場研究中，由于地震發震的快速性特點，單個地震的P、T軸分布往往不能作為真實發震區域的應力場主軸分布方位。數值模擬分析表明，當地震發生的數量足夠多時，大量震源機制解P、T軸的平均取向可代表該地區的主應力軸方位。圖2為本研究中相應區域的P、T軸平面分布圖，用以可視化日本至中國東北震源機制解的分區情況及檢驗應力反演結果的可靠性[11-12]。

圓形網格表示赤平投影結果，其中紅點表示P軸分布，藍點表示T軸分布圖2 震源機制解分區及各區P、T軸平面投影分布Fig.2 The distribution of focal mechanism solution and the plan projection of P axis and T axis of each region

2.2 各分區反演結果分析

2.2.1 深度0～30 km

矩陣[1,0,5]區經度為125°～130°，緯度為35°～40°，深度為0～30 km，本區反演共使用28個震源機制解，反演結果見圖2～3和表3～4。反演結果表明，σ1軸的最優方位角為-11.20°，角度變化范圍為-21.70°～0.10°，該區域主壓應力軸為NNW向，與菲律賓板塊擠壓歐亞板塊運動方向一致。傾伏角最優解為14.80°，角度變化范圍為-44.60°～26.40°，傾伏角較小。σ2軸方位角為-124.00°，為SSW向，傾伏角最優解為55.70°。σ3軸方位角最優解為87.60°，角度變化范圍為53.50°～248.40°，為NEE-NE向，張應力軸近水平，角度變化范圍較大，相應的最優解不確定度系數也較大。該區R值最優解為0.75，表明該區壓應力軸較穩定，中間應力軸與張應力軸性質趨于一致，反映該區壓應力軸表現明顯。

矩陣[2,0,5]區經度為130°～135°，緯度為35°～40°，深度為0～30 km，本區共使用240個震源機制解進行反演，反演結果見圖2～3和表3～4。反演結果表明，σ1軸方位角最優解為21.20°，角度變化范圍為17.00°～25.40°，為NNE向，主壓應力軸與矩陣[1,0,5]區相比存在明顯偏轉。σ3軸方位角最優解為118.00°，角度變化范圍為107.60°～116.00°，為SEE向。主張應力軸和主壓應力軸傾伏角近水平，R值最優解為0.45，表明該區主壓和主張應力軸均表現明顯。

緯度網格為Y軸，經度網格為X軸，深度網格為Z軸，紅色柱棒為σ1軸，綠色棒柱為σ2軸，藍色棒柱為σ3軸，紅色線段為Benioff帶走勢曲線圖3 各區網格矩陣構造應力場反演結果示意圖Fig.3 Diagram of inversion results of tectonic stress field of grid matrix in each region

矩陣[3,0,5]區經度為135°～140°，緯度為35°～40°，深度為0～30 km，該區反演共使用1 083個震源機制解， 反演結果見圖2～3和表3～4。反演結果表明，σ1軸最優方位角為25.30°，角度變化范圍為24.10°～28.50°，主壓應力軸仍為NNE向，表明在相同條件下，隨著經度的增加，主壓應力軸角度發生順時針偏轉。σ3軸方位角最優解為175.20°，角度變化范圍為133.40°～264.50°，為SSE-SE向，而該處σ3軸傾伏角最優解為82.80°，近于垂直。R值最優解為0.83，表明主張應力軸水平向出現硬性阻礙物，導致傾伏角和R值變大。

矩陣[3,1,5]區經度為135°～140°，緯度為40°～45°，深度為0～30 km，本區反演共使用137個震源機制解，反演結果見圖2～3和表3～4。反演結果表明，σ1軸方位角最優解為6.90°，角度變化范圍為-1.60°～12.00°，方位角為NNE-NE向，傾伏角最優解為6.20°，顯示出與[3,0,5]區相同的性質，傾伏角近水平。與[3,0,5]區的方位角相比，σ1軸隨緯度的增加，主壓應力軸趨于正NS向。σ3軸方位角最優解為-177.40°，角度變化范圍為-257.00°～-89.40°，方位角為NNW向，傾伏角最優解為83.80°，變化范圍為62.70°～88.20°。R值最優解為0.76，與同經度的[3,0,5]區相比，傾伏角表現結果一致，這進一步說明該區的應力場受某一區域塊體控制。

2.2.2 深度30～60 km

該深度僅矩陣[3,0,4]區滿足相關反演要求，其經度為135°～140°，緯度為35°～40°，深度為30～60 km，反演共使用276個震源機制解，反演結果見圖2～3和表3～4。反演結果表明，σ1軸最優方位角為62.60°，角度變化范圍為56.90°～70.60°，為NEE向，最優傾伏角較小，為23.30°。σ3軸方位角為171.70°，角度變化范圍為153.40°～218.00°，方位角變化較大，為SSE-SE向；傾伏角最優解為37.10°，變化范圍為-1.30°～66.00°。R值最優解為0.87，表明主壓應力軸方向作用明顯。

表3 各分區構造應力反演結果

表4 相關分區反演結果

2.2.3 深度60～100 km

該深度僅矩陣[3,0,3]區滿足反演所需的震源機制解個數要求，其經度為135°～140°，緯度為35°～40°，深度為60～100 km，該區反演共使用111個震源機制解，反演結果見圖2～3和表3～4。反演結果顯示，最大主壓應力軸方位角為26.70°，角度變化范圍為-8.30°～44.80°，為NNE-NE向，與相同經緯度、不同深度的[3,0,5]區和[3,0,4]區主壓應力軸基本保持一致。最小主應力軸最優方位角為122.10°，角度變化范圍為88.00°～285.90°，為SEE-SE向，角度變化范圍較大。σ3軸最優傾伏角為55.80°，角度變化范圍為-29.50°～79.00°。與[3,0,5]區和[3,0,4]區相比，主張應力軸傾伏角逐漸變小。R值最優解為0.78，表明該區仍被主壓應力軸所控制。

2.2.4 深度100～300 km

矩陣[3,0,2]區經度為135°～140°，緯度為35°～40°，深度為100～300 km，該區反演共使用169個震源機制解，反演結果見圖2～3和表3～4。反演結果表明，σ1軸最優方位角為-116.90°，角度變化范圍為-175.50°～48.30°，為SWW-SW向，與矩陣[3,0,5]、[3,0,4]和[3,0,3]區相比，偏轉角度逆時針陡增；σ1軸最優傾伏角為46.60°，角度變化范圍為-34.50°～87.40°，與同經緯度的[3,0,5]、[3,0,4]和[3,0,3]區相比，傾伏角明顯增大。σ3軸最優方位角為131.10°，角度變化范圍為78.90°～169.40°，最大主張應力軸以SSE向為主，最優傾伏角為19.50°，與同經緯度的[3,0,5]、[3,0,4]和[3,0,3]區相比，傾伏角明顯變緩，說明到達該深度時障礙物的擠壓作用明顯下降。該區R值最優解為0.32，表明該區主張應力作用較明顯。

矩陣[3,1,2]區經度為135°～140°，緯度為40°～45°，深度為100～300 km。該區σ1軸最優方位角為-153.80°，最大主壓應力軸為SSW向，最優傾伏角為46.80°，與同深度的[3,0,2]區基本一致，說明同深度主壓應力軸橫向壓應力保持穩定狀態。σ3軸最優方位角為-61.30°，最大主張應力軸為NWW向，傾伏角最優解為2.30°，與[3,0,2]區相比，角度偏轉較大，向相反方向發生張性破裂。R值最優解為0.89，表明該區以壓性應力狀態為主，與[3,0,2]區相反。

2.2.5 深度300～500 km

矩陣[3,0,1]區經度為135°～140°，緯度為35°～40°，深度為300～500 km，該區應力反演共使用75個震源機制解，反演結果見圖2～3和表3～4。反演結果表明，σ1軸最優方位角為167.80°，變化區間為162.00°～191.60°，主壓應力軸為SSE-SE向。與同經緯度、不同深度的[3,0,2]區存在明顯差異，說明在該區約300 km處存在使主壓應力軸發生偏轉的巨大作用力。σ1軸最優傾伏角為37.40°，傾伏角與[3,0,2]區相比明顯變緩，進一步說明該區分界處存在巨大作用力。σ3軸最優方位角為54.40°，變化范圍為38.80°～77.60°，最大主張應力軸為NEE向，傾伏角最優解為27.40°，角度變化范圍為5.80°～39.80°。該區張應力軸與同經緯度、不同深度的[3,0,2]區相比，NEE方向傾伏角變緩，表明阻礙張應力軸的障礙物逐漸減少。

矩陣[3,1,1]區經度為135°～140°，緯度為40°～45°，深度為300～500 km，該區反演共使用23個震源機制解。反演應力結果表明，壓應力軸與[3,0,1]區基本相同，說明該區的壓應力軸較穩定，張應力軸與[3,0,1]區正好相反，說明該深度的2個區域的張性破裂受到2股應力影響。

矩陣[2,0,1]區經度為130°～135°，緯度為35°～40°，深度為300～500 km，該區反演共使用30個震源機制解。反演結果顯示，σ1軸與鄰近矩陣區的軸方向呈現出一致的趨勢性擠壓，σ3軸也表現為趨勢性拉張。

2.2.6 深度500～700 km

該深度僅矩陣[2,1,0]區滿足要求，經度為130°～135°，緯度為40°～45°，深度為500～700 km，該區反演共使用20個震源機制解，反演結果見圖2～3和表3～4。反演結果表明，σ1軸方位角最優解為173.20°，角度變化范圍為153.80°～210.20°，傾伏角最優解為58.90°，最大主壓應力軸為SSE向。與300～500 km區域相比，方位角保持一定趨同性，傾伏角變大。σ3軸方位角最優解為-48.00°，變化范圍為-131.40°～56.90°，傾伏角最優解為24.40°，最大主張應力軸為NW向。與300～500 km區域相比，張應力軸向偏西方向偏轉，傾伏角趨于平緩。

3 結 語

本文采用3D構造應力場反演，通過空間分區的方法，將收集到的2 291個震源機制解分配到相應區域中。分析結果表明：

1)從日本東海岸到中國東北區域，構造應力場呈現明顯的分層現象，地殼表層應力場與深源應力場具有明顯差別，淺層構造應力場的主張和主壓應力隨經緯度的改變出現一定程度的旋轉變化，各區的壓應力軸受不同板塊控制。岡田義光[13]研究認為，日本東北部屬于北美板塊，且日本東部邊緣正在形成新的板塊邊界。由此可以看出，地殼淺層中各區主壓應力與張應力不僅受到歐亞板塊、菲律賓板塊和太平洋板塊的相互擠壓碰撞，同時也受北美板塊的斜插擠壓作用，地殼應力結構復雜。矩陣[1,0,5]區的主壓應力軸為NNW向，說明該區域的壓應力主要來自于菲律賓板塊，太平洋板塊向西擠壓為輔。[2,0,5]區與[3,0,5]區的主壓應力軸為NNE向，與太平洋板塊擠壓歐亞板塊的方向存在偏差，進一步說明該區處于板塊交接處，受到北美板塊強大的斜插作用。[3,1,5]區壓應力軸方位角為6.90°，與[2,0,5]區和[3,0,5]區相比，壓應力軸呈近NS向，說明在經度相同的0～30 km深度，緯度越高，壓應力軸受到的偏轉作用力越強，北美板塊的擠壓斜插作用越劇烈。

2)由表2和圖3可知，隨著深度的增加，壓應力軸及張應力軸在約300 km深處出現明顯變化，最大主壓應力軸由水平方向突然向下俯沖，隨后俯沖速度以一定規律變緩，說明該區域在300 km以下已到達Benioff俯沖帶，且俯沖帶一直延伸到700 km左右。

3)由表2和圖3可以看出，最大主張應力軸在北美板塊與歐亞板塊的交接地帶呈現出一定規律的俯沖現象，俯沖由淺部一直持續到300～500 km左右并逐漸變緩消亡。可以看出，隨著時間的推移，北美板塊在一定程度上繼續向NE方向擠壓歐亞板塊，并持續發生地震。

4)由三維構造應力場圖中主壓應力軸可以看出，隨著深度的增加，Benioff帶向下凹陷俯沖。圖4為旋轉后的三維立體圖，從圖中可以看出，俯沖過程一直延續到中國東北區域，該區域的深震與Benioff帶的俯沖作用具有重要聯系，Benioff帶沿著日本海溝一線斜插入中國東北區域約500～700 km深處。

紅色散點為σ1軸分布，綠色散點為σ2軸分布，藍色散點為σ3軸分布，紅色線段為Benioff帶走勢曲線圖4 3D網格矩陣旋轉后各區散點構造應力場反演結果示意圖Fig.4 Diagram of inversion results of regional scattered tectonic stress field after 3D grid matrix rotation

本次研究以空間三維的角度進行構造應力場反演，可更清晰有效地對某一地區的地塊運動進行判斷分析，但變化無處不在。本研究尚未涉及以時間軸為變化的空間四維應力場反演，后續將在該領域繼續探索。