孟國杰， 蘇小寧， 徐婉楨， 任金衛(wèi)， 楊永林, Nikolay V. Shestakov

1 中國地震局地震預(yù)測重點實驗室(中國地震局地震預(yù)測研究所)， 北京 100036 2 四川省地震局測繪工程院， 四川雅安 625000 3 Institute of Applied Mathematics, Far Eastern Branch of Russian Academy of Sciences, Vladivostok 690091, Russia 4 Far Eastern Federal University, Vladivostok 690068, Russia

?

基于GPS觀測研究2010年青海玉樹MS7.1地震震后地殼形變特征及其機(jī)制

孟國杰1， 蘇小寧1， 徐婉楨1， 任金衛(wèi)1， 楊永林2, Nikolay V. Shestakov3,4

1 中國地震局地震預(yù)測重點實驗室(中國地震局地震預(yù)測研究所)， 北京 100036 2 四川省地震局測繪工程院， 四川雅安 625000 3 Institute of Applied Mathematics, Far Eastern Branch of Russian Academy of Sciences, Vladivostok 690091, Russia 4 Far Eastern Federal University, Vladivostok 690068, Russia

2010年4月14日青海玉樹MS7.1地震發(fā)生在青藏高原東南部甘孜—玉樹地震帶，在震后7～10天內(nèi)，我們快速建立了由15個GPS測站組成的跨地震破裂帶觀測剖面，包括1個連續(xù)站，3個半連續(xù)站和11個流動站，對所有站進(jìn)行了240多天的觀測，獲取了該次地震的震后形變時空特征.采用歐拉矢量和位錯模型解算了背景速度場，并從GPS觀測的形變場中扣除該分量.采用分層黏彈性位錯模型計算余震引起的地表形變，結(jié)果表明余震對部分測站的位移造成不可忽視的影響.采用對數(shù)模型擬合位移時間序列，表明特征衰減時間為6.7±1.2天.利用最速下降法反演震后余滑時空分布，反演結(jié)果表明震后斷層活動以左旋滑動為主，斷層南盤具有少量的抬升.在空間分布上，余滑主要位于同震破裂區(qū)的兩側(cè)，西北側(cè)的余滑幾乎達(dá)到地表，而東南區(qū)的余滑基本在同震破裂區(qū)的下方，余滑最大的區(qū)域位于結(jié)古鎮(zhèn)東南下方10～20 km的深度范圍.隨著震后離逝時間的增加，2個余滑區(qū)在空間上保持不變，余滑區(qū)的面積逐漸擴(kuò)大.余滑的矩釋放為(1.5～5.1)×1018Nm ，相當(dāng)于1個MW6.1～6.4地震釋放的能量.分層巖石圈黏彈性模型計算的地殼孔隙彈性反彈形變與地表觀測值相差較大，不能解釋觀測到的震后變形.采用麥克斯維爾流變體模型計算下地殼和上地幔松弛引起的地表形變，顯示出其對地表形變的貢獻(xiàn)較小.GPS觀測得到的震后形變所具有的快速衰減特征，以及余滑模型能夠較好地擬合GPS地表形變，表明2010年玉樹MS7.1地震后早期階段的地殼形變主要是由余滑機(jī)制決定的.

玉樹7.1級地震； 震后形變； 構(gòu)造活動； 余震活動； 震后余滑

1 引言

2010年4月14日7時49分(當(dāng)?shù)貢r間)青海省玉樹藏族自治州發(fā)生MS7.1地震，震中位于(33.11°N，96.59°E)，震源深度14 km(中國地震臺網(wǎng)統(tǒng)一正式地震目錄，http:∥10.5.202.22).玉樹地震發(fā)生在甘孜—玉樹斷裂帶上，該斷裂帶是一條沿北西—南東方向延伸、以左旋走滑為主的斷裂帶，是巴顏喀拉地塊和羌塘地塊的構(gòu)造邊界，總長約500 km(聞學(xué)澤等，2003).此次地震在地表引起一條長約50～65 km的破裂帶，以走滑為主，具有50 m的寬度(Chen et al.，2010; Li et al., 2012).發(fā)震斷層面走向119°，傾角83°，滑動角-2°，為高角度左旋走滑破裂(張勇等，2010).自主震發(fā)生時刻至12月31日共發(fā)生余震2千多次，最大震級達(dá)到6.4級(王未來等，2012).

利用GPS等大地測量技術(shù)可以觀測大地震之后的地表形變，獲取震后形變特征(任金衛(wèi)等，2005；丁開華等，2013)，進(jìn)而分析震后形變機(jī)制、研究斷層摩擦性質(zhì)和巖石圈黏彈性結(jié)構(gòu)(Wen et al., 2012；Gonzalez-ortega et al., 2014).GPS觀測的震后形變可能包括有多種因素的貢獻(xiàn)，如區(qū)域構(gòu)造運動、余震活動、孔隙彈性反彈、余滑、巖石圈應(yīng)力松弛等(Wang et al.,2010).前人研究結(jié)果顯示震后余滑和巖石圈應(yīng)力松弛是震后形變的主要因素(Marone et al. 1991; Shen et al.,1992; Pollitz et al.,1997)，在震后早期階段震后余滑占主導(dǎo)地位，在震后的中期至后期階段，巖石圈應(yīng)力松弛占主導(dǎo)地位(Wang et al., 2009).

玉樹地震發(fā)生后，中國地震局地震預(yù)測研究所聯(lián)合四川省地震局測繪工程院立即成立應(yīng)急觀測小組赴震區(qū)開展GPS建站和觀測工作，在震后7～10天內(nèi)新建9個GPS測站，利用6個已有的“中國地殼運動觀測網(wǎng)絡(luò)”測站，組成一條跨玉樹發(fā)震斷裂的觀測剖面(圖1)，剖面長度360 km左右，其中包括1個連續(xù)觀測站(QHYS)、3個半連續(xù)觀測站(J392，YS01, YS04)和11個流動觀測站，對該剖面進(jìn)行了240多天的觀測.本文通過對觀測數(shù)據(jù)進(jìn)行高精度解析得到相對ITRF2005(Altamimi et al., 2007)的測站坐標(biāo)時間序列，利用構(gòu)造塊體歐拉矢量模型和孕震層下方脆韌轉(zhuǎn)換帶蠕動模型去除構(gòu)造活動引起的背景運動.基于雙差定位余震目錄，計算余震引起的形變，并對其進(jìn)行改正.采用對數(shù)模型對位移時間序列進(jìn)行擬合，給出了震后地表形變時空衰減特征.采用最速下降法反演了震后余滑分布，并對震后形變機(jī)制進(jìn)行了分析和討論.

圖1 (a) 研究區(qū)域位置; (b) 青海玉樹區(qū)域地震構(gòu)造(a) 圖中黑色方框代表研究區(qū)域，黑色實線為青藏高原活動地塊邊界.藍(lán)色五角星為玉樹地震震中位置; (b) 圖中“沙灘球”表示玉樹地震震源機(jī)制解，紅色圓為余震分布，黑色實線表示斷層.Fig.1 (a) Location of study area; (b) Seismo-tectonic of Yushu area, Qinghai Province

2 GPS觀測與數(shù)據(jù)處理

在2010年4月20日完成GPS觀測儀器的架設(shè)后，青海玉樹站(QHYS)作為唯一的連續(xù)站一直保持連續(xù)觀測.3個半連續(xù)站和和11個流動觀測站自2010年4月20日至11月25日共觀測5期(表1)，第一期包括建站和觀測，并持續(xù)11天，其后4期觀測每期都為6天.基于震后地表形變衰減快的特點確定相鄰2期觀測間隔時間，越到后期，間隔時間越長.5期觀測均采用Trimble 5700、Trimble NetRS GPS接收機(jī)和TRM29659.00扼徑圈天線，在各期觀測中保證在同一站點采用同一臺接收機(jī)和天線，以減少因接收機(jī)和天線不同可能導(dǎo)致的系統(tǒng)誤差.

表1 玉樹地震震后地殼形變GPS觀測方案

采用GAMIT/GLOBK軟件(Herring et al., 2008)進(jìn)行GPS數(shù)據(jù)處理，GPS數(shù)據(jù)包括本文15個剖面站、“網(wǎng)絡(luò)工程”25個基準(zhǔn)站以及中國大陸周邊13個IGS連續(xù)站的數(shù)據(jù).對衛(wèi)星軌道進(jìn)行一定的松弛，涉及的地球重力場、固體潮和極潮模型遵循IERS2003規(guī)范.采用FES2004模型改正海潮引起的測站形變.選擇衛(wèi)星截止高度角大于10°的數(shù)據(jù)用于數(shù)據(jù)編輯和參數(shù)估計，觀測數(shù)據(jù)的定權(quán)采用衛(wèi)星截止高度角的函數(shù)，采用LC組合觀測量以消除電離層一階項影響.對流層天頂延遲改正采用GPT模型，對每個測站每1個小時估計1個天頂延遲參數(shù)，映射函數(shù)采用GMF模型，對每個測站的EW向和NS向各附加1個大氣水平梯度參數(shù).把解算的單日松弛解和全球 IGS 核心站單日松弛解合并，得到包含剖面站、“網(wǎng)絡(luò)工程”基準(zhǔn)站和IGS 核心站的單日松弛解.選擇較均勻地分布于全球、在ITRF2005全球參考框架(Altamimi et al., 2007)中水平速度分量中誤差較小的IGS站作為框架站，并將IGS核心站的坐標(biāo)和速度約束于ITRF2005框架，解算GPS剖面站在ITRF2005中的坐標(biāo)，進(jìn)而通過坐標(biāo)差分得到各測站相對于發(fā)震時刻的位移時間序列.

3 構(gòu)造運動引起的地殼背景形變及其改正

經(jīng)過上述GPS數(shù)據(jù)處理得到的相對于ITRF2005的位移時間序列包含構(gòu)造的背景運動和震后形變，為了提取出震后形變，需從中扣除構(gòu)造運動引起的位移分量.構(gòu)造運動引起的背景運動可由地震前多年的GPS觀測得到，由于本文的GPS測站有9個是震后7～ 10天內(nèi)建立的新站，6個先存的“網(wǎng)絡(luò)工程”區(qū)域站在震前只有一期觀測，因此，這15個站點都缺乏背景速度.

本文收集了青藏高原地區(qū)“中國地殼運動觀測網(wǎng)絡(luò)”自1999年至2010年玉樹MS7.1地震前所有GPS原始觀測數(shù)據(jù)，并進(jìn)行高精度數(shù)據(jù)處理，獲取了相對于ITRF2005(Altamimi et al., 2007)的速度場.GPS站的背景運動速度包括2部分： (1) 彈性塊體的長期運動速度； (2) 孕震層下方脆韌轉(zhuǎn)換帶蠕動引起的形變.前者可由歐拉矢量模型計算(Argus et al., 2010; Altamimi et al., 2012)，后者可由位錯模型計算(Savage & Burford，1973).考慮到甘孜—玉樹斷層具有較強(qiáng)的活動性(聞學(xué)澤等，2003)，跨該斷裂的研究區(qū)不可以作為一個整體計算其歐拉矢量，本文分別計算斷層兩側(cè)巴顏喀拉塊體和羌塘塊體相對ITRF2005的歐拉矢量，然后用之分別計算發(fā)震斷層?xùn)|北側(cè)和西南側(cè)GPS站的速度，計算的EW速度在41.7～47.9 mm·a-1之間，NS速度在-8.98 ～-2.11 mm·a-1之間.采用線性關(guān)系計算每一個觀測時間的位移，并從各GPS測站相對于ITRF2005的位移時間序列中扣除該分量.

對于孕震層下方脆韌轉(zhuǎn)換帶蠕動引起的地表形變，可以采用彈性半空間的彈性模型(Savage & Burford, 1973)，也可以采用水平分層的黏彈性位錯模型(Wang et al., 2006)進(jìn)行計算，本文采用后者進(jìn)行計算.鑒于甘孜—玉樹斷裂帶上地震震源深度分布主要位于5～20 km(王未來等，2012)，因此假定地表至地下20 km為完全閉鎖.為了避免斷層端部效應(yīng)的影響，假定在20 km以下存在一條長度為400 km，深度方向?qū)挾葹?00 km的位錯，走向為119°，上邊界的中點位于玉樹地震震源.采用Ryder等(2014)的巖石圈分層模型，參考玉樹地震的震源機(jī)制(張勇等，2010)，設(shè)位錯的傾角為83°，滑動方式為左旋錯動.參考Lin等(2011)研究結(jié)果，設(shè)脆韌轉(zhuǎn)換帶蠕動速率為8 mm·a-1.計算結(jié)果(表2)顯示，影響最大的站為YS07，其EW和NS速度分別為0.70～0.38 mm·a-1，總速度為0.80 mm·a-1.由此可知，孕震層下方脆韌轉(zhuǎn)換帶蠕動引起的地表形變較小，因此在后續(xù)的分析中不予考慮.

表2 脆韌轉(zhuǎn)換帶蠕動引起的GPS站背景速度

4 余震引起的地表形變

本文挑選14個4級以上的余震，計算其在地表引起的位移.假定這些余震的震源機(jī)制與主震的相同，余震破裂面的尺度及滑動量由尺度定理(Wells 和Coppersmith, 1994)進(jìn)行計算.采用Ryder等(2014)給出的巖石圈黏彈性分層模型，上層為彈性層，其厚度為15 km，下層為Maxwell黏彈性層，其黏滯度為6×1019Pa·s.利用矩形位錯模型計算余震在地表15個GPS測站引起的位移(圖2)，應(yīng)用PSGRN/PCMP程序進(jìn)行計算(Wang et al., 2006).

圖2表明主震后1.5 h發(fā)生的MS6.4余震和第46天發(fā)生的MS5.9余震影響較大，其余震級較小的余震引起的位移量很小.余震活動在一些測站引起明顯的地表位移，余震導(dǎo)致YS05站EW位移達(dá)到約3 mm，NS位移達(dá)到約7 mm.由于余震對GPS觀測的震后形變影響不可忽略，需要從GPS觀測的位移時間序列中扣除其分量.需要說明的是，由于YS03、J390、J396、J397、J399測站NS位移值與其他測站接近，YS01的覘標(biāo)曾遭受一定的破壞，為使圖形清晰，余震在這些站引起的位移未在圖2中畫出.

圖2 M≥4.0余震及其在GPS測站引起的地表形變橫軸坐標(biāo)原點為主震發(fā)生時刻，折線中的階變表示余震引起的同震地表位移.Fig.2 M≥4 aftershocks and induced surface deformation at GPS sitesOrigin of horizontal axis is the occurrence time of mainshock. Steps in the broken lines indicate displacements due to the aftershocks.

圖3 玉樹斷層北側(cè)GPS測站震后位移對數(shù)模型擬合結(jié)果(實心圓圈為觀測位移，藍(lán)線為模型曲線)Fig.3 Logarithmic fitting of postseismic displacementsfor GPS sites north of Yushu fault (Solid circles indicate observational displacements, and blue curves are modeled displacements)

圖4 玉樹斷層南側(cè)GPS測站震后位移對數(shù)模型擬合結(jié)果(實心圓圈為觀測位移，藍(lán)線為模型結(jié)果)Fig.4 Logarithmic fitting of postseismic displacements for GPS sites south of Yushu fault(Solid circles indicate observational displacements, and blue curves are modeled displacements)

5 GPS觀測的震后形變與函數(shù)擬合

對經(jīng)過構(gòu)造形變、余震形變改正的位移時間序列進(jìn)行函數(shù)擬合，以分析震后形變的連續(xù)變形過程.可以采用不同的函數(shù)模型，如指數(shù)衰減模型(Shen et al., 1994)、對數(shù)衰減模型(Savage and Svarc, 2009)和速率增強(qiáng)的摩擦律(Barbot et al.，2009)對震后形變進(jìn)行擬合.本文采用如下對數(shù)模型擬合震后位移：

(1)

其中，a為震后位移時間序列初始常數(shù)，b為震后位移衰減振幅，teq為地震發(fā)生的時間，τ為衰減時間.式中a、b和τ為待估參數(shù).對τ全局進(jìn)行搜索，把非線性方程轉(zhuǎn)化為線性方程，采用奇異值分解法對法方程進(jìn)行求逆，并且約束EW和NS分量具有相同的衰減時間.由于各站總體觀測時間長度不盡相同，故選擇其中7個總體觀測時間較長，觀測質(zhì)量較好的QHYS、J396、YS08、J399、YS06、YS05和YS07測站，對其位移時間序列進(jìn)行擬合，解算出的特征衰減時間為6.7±0.2天，然后將其代入公式(1)，分別解算其參數(shù)a和b.對于余下的7個測站，固定其特征衰減時間，分別計算其常數(shù)a和振幅b，計算結(jié)果見表3.

圖3和圖4分別給出了玉樹斷層北側(cè)和南側(cè)GPS測站震后變形時間序列及對數(shù)模型擬合結(jié)果.總體上顯示出測站越靠近斷層，其震后變形越大，隨著測站遠(yuǎn)離斷層，其變形幅度逐漸減少.由于測站相對稀少，無法測試變形是否相對斷層破裂帶具有對稱性，但斷層兩側(cè)測站的位移大致是相等的.表4給出了震后250天各測站累計的震后位移，該表顯示出距斷層最近的QHYS站具有最大的形變，震后250天累積位移達(dá)到17.57 mm.其次是位于破裂帶東北側(cè)的YS04和J397站，二者都表現(xiàn)為西向運動為主，兼有有少量南向運動，位移幅度分別達(dá)到15.82 mm和14.50 mm.距離斷層最遠(yuǎn)的測站J400和J390站分別位于斷層的最北側(cè)和最南側(cè)，變形幅度分別為2.64 mm和2.48 mm，這兩個測站距離斷層破裂帶的距離都超過175 km.位于斷層以北的測站具有少量的南向運動，而位于斷層以南的測站具有少量的北向運動.總體上顯示出主震發(fā)生后玉樹斷層具有以左旋為主，兼有少量擠壓的特征.

主震之后的前50天，各測站的形變顯示出快速的衰減過程，之后衰減變慢，震后前50天的位移量與震后250天的位移量的比值都超過50%.就玉樹連續(xù)站而言，震后50～100、100～150、150～200天、200～250天的位移量與震后250天內(nèi)位移的比值分別為58.61%，17.33%，10.58%,5.93%，顯示出比較完整的衰減過程.

表3 震后位移時間序列對數(shù)函數(shù)擬合結(jié)果

表4 由對數(shù)模型計算的震后250天累積位移

6 震后余滑反演

前人基于InSAR、GPS、地震波形數(shù)據(jù)反演了玉樹地震的同震破裂模型(Li et al., 2011; Tobita et al., 2011; Jiang et al., 2013; Meng et al., 2013; Qu et al., 2013，李奇, 2015).本文參考李奇(2015)聯(lián)合GPS和InSAR數(shù)據(jù)的反演結(jié)果，把斷層面離散化為規(guī)則的子斷層，每個子斷層的尺度為2 km×2 km.仍采用分層的黏彈性巖石圈模型(Ryder, et al., 2014).利用最速下降法(Cervelli et al.，2001)反演各子斷層的滑移矢量.采用拉普拉斯算子平滑相鄰的子斷層位錯量，該算子的平滑作用由平滑因子控制，平滑因子α的選擇由模型粗糙度相對數(shù)據(jù)擬合殘差的L曲線確定(Segall and Harris 1987； D′Agostino et al., 2012).

以本文第5節(jié)給出的對數(shù)模型計算震后位移，并依此為約束數(shù)據(jù)反演震后余滑，采用位移誤差平方的倒數(shù)作為其權(quán)重，對于土層站點的位移數(shù)據(jù)，作進(jìn)一步降權(quán)，使其權(quán)重降為基巖點的60%.通過改變數(shù)據(jù)自由度進(jìn)行多次反演.首先約束斷層滑動為左旋走滑，然后假定斷層為走滑兼逆沖，反演結(jié)果顯示出觀測數(shù)據(jù)主要可以由左旋走滑分量進(jìn)行解釋，對傾滑分量的方向不進(jìn)行約束，容許其沿著斷層發(fā)生變化.對走滑和傾滑分量施加平滑約束，并采用L曲線方法確定數(shù)據(jù)的擬合質(zhì)量和模型的粗糙度.在反演過程中調(diào)整平滑因子，得到L曲線如圖5所示.平滑因子α取為0.1，反演得到震后250天余滑的空間分布(圖6).由觀測數(shù)據(jù)與模型位移的對比(圖7)可知，所有測站的觀測數(shù)據(jù)與模型計算的位移具有較好的一致性，大部分位移的殘差小于2.0 mm，數(shù)據(jù)和模型的相關(guān)性達(dá)到97.8%.

圖5 歸一化數(shù)據(jù)擬合殘差與模型粗糙度歸一化值關(guān)系模型粗糙度表示余滑在空間上的梯度，歸一化殘差表示擬合殘差平均值的相對大小.Fig.5 Relationship between normalized model roughness and normalized residualsModel roughness indicates the roughness of model slip; Normalized fitting residual indicates the average relative residuals.

余滑的空間分布顯示出，在震后250天斷層的余滑角2～10°，表明斷層活動以為左旋滑動為主，斷層南盤具有少量的抬升.與同震滑動分布模型對比可知，余滑與同震破裂在空間分布上具有一定的互補(bǔ)性，余滑主要位于同震破裂區(qū)的兩側(cè)，西北側(cè)的余滑達(dá)到地表，最大滑動量為0.19 m.東南端的余滑區(qū)位于同震破裂區(qū)的下方，深度在10～20 km范圍，最大余滑量達(dá)到0.23 m.屈春燕等(2013)基于ENVISAT/ASAR雷達(dá)數(shù)據(jù)給出了玉樹震后的形變場，發(fā)現(xiàn)在結(jié)古鎮(zhèn)西北方向主震的微觀震中區(qū)域存在較大范圍的形變，并認(rèn)為是震后余滑引起的.本文反演得到的西北端的余滑區(qū)與其得到的地表形變存在一致性.

圖6 余滑與同震位錯對比紅色實心圓圈為余震，黃色五角星為主震震源，黑色箭頭為余滑矢量.(a) 同震位錯模型(李奇，2015)； (b) 主震后250天余滑分布.Fig.6 Comparison of after slip and coseismic slipRed solid circles indicate aftershocks, yellow asterisk is focus of mainshock, and black arrows are afterslip vector. (a) Coseismic slip distribution; (b) After slip during 250 days after the mainshock.

對斷層面上余滑量進(jìn)行積分，并乘以地殼剪切模量可以計算出震后余滑的距釋放.對于破裂的地殼巖石，假定其剪切模量為10～33 GPa(Hamiel et al., 2006；Cochran et al., 2009).在震后250天內(nèi)，震后余滑的矩釋放為(1.5～5.1)×1018Nm，余震釋放的矩為2.2×1017Nm，余滑矩釋放與余震距釋放的比值7.0～22.7，顯示出玉樹地震震后形變主要是非震形變引起的. 利用矩震級計算公式(Kanamori, 1977)：

(2)

式(2)中，M0為余滑釋放的距，MW為矩震級，可計算出玉樹地震余滑相當(dāng)于釋放一次MW6.1～6.4地震的能量.

利用對數(shù)擬合模型計算震后50、100、150天和200天4個時間段的累計位移，并以此為約束反演震后余滑的時空演化(圖8)，發(fā)現(xiàn)西北側(cè)和東南側(cè)2個余滑區(qū)在空間位置上保持不變，隨著離逝時間增加，余滑的區(qū)域逐漸擴(kuò)大.

7 討論

已有研究表明中強(qiáng)以上地震之后的早期階段，地殼形變速率高于中期和后期階段(Ozawa et al., 2012; Wang et al., 2012；Sun et al., 2014)，形變速率的快速降低是震后早期形變的基本特征(Hu and Wang, 2012; Sun et al., 2014).由于巖石圈對同震應(yīng)力變化的響應(yīng)，地震后產(chǎn)生瞬態(tài)形變，其物理成因主要包括如下幾種：震后余滑(Marone et al., 1991; Perfettini et al., 2004; Barbot et al., 2009)，液體飽和的地殼孔隙彈性反彈(Peltzer et al., 1998，Jonsson et al., 2003; Fialko, 2004), 下地殼與上地幔中的黏彈性松弛，斷層區(qū)體積變化(Massonnet et al., 1996; Jacobs et al., 2002; Fielding et al., 2009).GPS接收機(jī)記錄的震后形變是多種因素的綜合影響結(jié)果，目前利用震后形變數(shù)據(jù)完全分離各種因素的影響或評估其相對重要性仍然是十分困難的(Hu and Wang, 2012; Nykolaishen, 2015)，尤其是僅具有震后早期階段短期的形變時間序列(Wang et al., 2012).

圖7 主震后250天地表累積位移和余滑模型計算的位移對比(“沙灘球”為主震的震源機(jī)制)Fig.7 Cumulative displacements during 250 days following the mainshock and the displacement predicted by the afterslip model (The beachball indicates focal mechanism of the mainshock)

7.1 孔隙彈性反彈引起的地表形變

地殼孔隙彈性反彈是由震時孔隙壓力改變造成孔隙流體滲流引起的，孔隙流體反彈時間尺度依賴于巖石的有效滲透率，通常持續(xù)數(shù)天到數(shù)年(Peltzer et al.,1998；Hu et al., 2014)，其引起的地表形變難以由解析方法得到，其時間衰減是非線性的(Fialko, 2004).但可利用地殼的2個狀態(tài)估計孔隙彈性反彈引起的地表形變，第一個狀態(tài)代表地殼對地震的瞬時響應(yīng)，稱為未排空水(undrained)狀態(tài).第二個狀態(tài)是孔隙流體壓力引起的擾動達(dá)到平衡條件，短期的孔隙彈性流體流動已經(jīng)完成，達(dá)到平衡狀態(tài).這與排干狀態(tài)(drained)不同，排干狀態(tài)表示由于緩慢的加載過程和較高的可滲透率，孔隙流體壓力沒有變化(Hu et al.,2014).利用同震破裂模型計算這兩種狀態(tài)下形變的差異可以近似為時間相依的孔隙彈性反彈總的影響，代表孔隙彈性反彈形變的最大估計.液體飽和的地殼孔隙彈性反彈引起的地表形變主要集中在發(fā)震斷層帶的端部(Peltzer et al., 1998; Jonsson et al., 2003).

利用地震前后泊松比的變化，采用均勻彈性半空間或分層黏彈性位錯模型可模擬孔隙彈性反彈引起的累計震后形變(Segall, 2010).假定未排干流體的地殼泊松比υu= 0.27，平衡狀態(tài)的地殼泊松比υd=0.31(Riceand Cleary 1976；Peltzer et al., 1998).采用Ryder 等(2014)的巖石圈模型和李奇(2015)的同震滑動模型，以兩種泊松比分別計算地表位移場，并對其進(jìn)行差分得到平衡狀態(tài)的孔隙彈性響應(yīng)形變(圖9).圖9顯示出孔隙彈性模型給出的地表形變與GPS觀測的形變沒有一致性，無法用該模型解釋震后變形.考慮到玉樹地區(qū)為比較干旱，可能不處于飽和水的狀態(tài)，給定的孔隙彈性反彈模型計算的位移與GPS觀測的位移沒有相似性.因此，可以認(rèn)為地殼孔隙彈性反彈引起的地表位移貢獻(xiàn)量較小.

圖8 余滑時空演化紅色箭頭表示余滑矢量，五角星為主震的震源. (a) 震后50天； (b) 震后100天； (c) 震后150天； (d) 震后200天.Fig.8 Spatiotemporal afterslip distributionRed arrows indicate afterslip vectors，and the asterisk is focus of mainshock. (a), (b), (c), (d) are afterslips in the periods of 50, 100, 150 and 250 days following the main shock, respectively.

圖9 地殼孔隙彈性反彈引起的地表形變紅色箭頭表示GPS觀測的震后位移，黑色箭頭表示模型結(jié)果,藍(lán)線為地表破裂.Fig.9 Surface deformation induced by poroelastic reboundRed and black arrows indicate observed and modelled postseismic displacements, respectively.

圖10 下地殼與上地幔黏彈性松弛引起的地表形變Fig.10 Surficial deformation caused by viscoelastic relaxation of lower crust and upper mantle

7.2 下地殼與上地幔黏彈性松弛引起的地表形變

同震引起的應(yīng)力變化對黏彈性的下地殼與上地幔施加力的作用，因此巖石圈發(fā)生松弛效應(yīng)引起地表發(fā)生形變，并可持續(xù)多年(Calais et al., 2002；Wang et al., 2012； Ryder et al., 2014).采用Ryder等(2014)給出的巖石圈黏彈性模型和李奇(2015)同震滑動模型，計算玉樹地震后巖石圈黏彈性松弛引起的地表形變.計算結(jié)果(圖10)顯示出下地殼與上地幔黏彈性松弛對GPS測站引起的影響較小，在震后250天內(nèi)累計黏彈性形變最大的測站為YS03，其EW向形變量值為0.21 mm.

8 結(jié)論

2010年青海玉樹MS7.1地震發(fā)生后，我們迅速在震中區(qū)建立15個GPS測站組成的跨斷層觀測剖面，并進(jìn)行了240多天的觀測，獲取了該次地震的震后形變特征.在震后的前50天地殼形變量大，衰減迅速，所有測站在此期間的位移量都達(dá)到震后250天理論位移量的58%.余震對部分測站位移具有不可忽視的影響，地殼介質(zhì)孔隙彈性反彈模型不能解釋GPS觀測的震后變形，巖石圈粘彈性松弛引起的震后變形可以忽略不計.

余滑反演結(jié)果表明玉樹地震后發(fā)震斷層以左旋滑動為主，斷層南盤具有少量的抬升.余滑與同震破裂在空間分布上具有一定的互補(bǔ)性，余滑主要位于同震破裂區(qū)的兩端，西北端的余滑區(qū)幾乎達(dá)到地表，而東南區(qū)的余滑基本在同震破裂區(qū)的下方，深度在10～20 km范圍.隨著震后離逝時間的增加，兩個余滑區(qū)在空間上保持穩(wěn)定，面積逐漸擴(kuò)大.GPS觀測獲得的震后形變所具有的快速衰減特征，以及余滑模型能夠較好地擬合震后地表形變，表明2010年MS7.1玉樹地震后地表形變主要是由余滑機(jī)制決定的.

致謝 中國地震局地震預(yù)測研究所正研級高工孫漢榮設(shè)計了GPS建站方案，中國地震局地震預(yù)測研究所、四川省地震局測繪工程院和青海省地震局多位同仁參加了玉樹地震GPS野外建站與觀測工作，德國地球科學(xué)中心(GFZ)汪榮江在震后余滑和巖石圈應(yīng)力松弛模擬方面提供熱情幫助，石玉濤博士提供雙差定位余震目錄，兩位審稿者提出了寶貴的修改意見，作者在此一并致謝.

Altamimi Z, Collilieux X, Legrand J, et al. 2007. ITRF2005: A new release of the international terrestrial reference frame based on time series of station positions and Earth orientation parameters.J.Geophys.Res., 112(B9): B09401, doi: 10.1029/2007JB004949.Altamimi Z, Métivier L, Collilieux X. 2012. ITRF2008 plate motion model.J.Geophys.Res., 117(B7): B07402, doi: 10.1029/2011JB008930. Argus D F, Gordon R G, Heflin M B, et al. 2010. The angular velocities of the plates and the velocity of Earth′s centre from space geodesy.Geophys.J.Int., 180(3): 913-960.

Barbot S, Fialko Y, Bock Y. 2009. Postseismic deformation due to theMW6.0 2004 Parkfield earthquake: Stress-driven creep on a fault with spatially variable rate-and-state friction parameters.J.Geophys.Res., 114(B7): B07405, doi: 10.1029/2008JB005748.

Calais E, Vergnolle M, Déverchère J, et al. 2002. Are post-seismic effects of theM=8.4 Bolnay earthquake (1905 July 23) still influencing GPS velocities in the Mongolia-Baikal area?Geophys.J.Int., 149(1): 157-168.Cervelli P, Murray M H, Segall P, et al. 2001. Estimating source parameters from deformation data, with an application to the March 1997 earthquake swarm off the Izu Peninsula, Japan.J.Geophys.Res., 106(B6): 11217-11237.

Chen L C, Wang H, Ran Y K, et al. 2010. TheMS7.1 Yushu earthquake surface rupture and large historical earthquakes on the Garz-Yushu fault.Chin.Sci.Bull., 55(31): 3504-3509.

Cochran E S, Li Y G, Shearer P M, et al. 2009. Seismic and geodetic evidence for extensive, long-lived fault damage zones.Geology, 37(4): 315-318.

D′Agostino N, Cheloni D, Fornaro G, et al. 2012. Space-time distribution of afterslip following the 2009 L′Aquila earthquake.J.Geophys.Res., 117(B2): B02402, doi: 10.1029/2011JB008523. Ding K H, Xu C J, Wen Y M. 2013. Postseismic deformation associated with the 2008 Wenchuan earthquake by GPS data.GeomaticsandInformationScienceofWuhuaUniversity(in Chinese), 38(2): 131-135.

Fialko Y. 2004. Evidence of fluid-filled upper crust from observations of postseismic deformation due to the 1992MW7.3 Landers earthquake.J.Geophys.Res., 109(B8): B08401, doi: 10.1029/2004JB002985. Fielding E J, Lundgren P R, Bürgmann R, et al. 2009. Shallow fault-zone dilatancy recovery after the 2003 Bam earthquake in Iran.Nature, 458(7234): 64-68.

Gan W J, Zhang P Z, Shen Z K, et al. 2007. Present-day crustal motion within the Tibetan Plateau inferred from GPS measurements.J.Geophys.Res., 112(B8): B08416, doi: 10.1029/2005jb004120. Gonzalez-Ortega A, Fialko Y, Sandwell D, et al. 2014. El Mayor-Cucapah (MW7.2) earthquake: Early near-field postseismic deformation from InSAR and GPS observations.J.Geophys.Res., 119(2): 1482-1497, doi: 10.1002/2013JB010193.

Hamiel Y, Katz O, Lyakhovsky V, et al. 2006. Stable and unstable damage evolution in rocks with implications to fracturing of granite.Geophys.J.Int., 167(2): 1005-1016.

Herring T A, King R W, Floyd M A, et al. 2008. Introduction to GAMIT/GLOBK. Cambridge, Mass: Massachusetts Institute of Technology.

Hu Y, Bürgmann R, Freymueller J T, et al. 2014. Contributions of poroelastic rebound and a weak volcanic arc to the postseismic deformation of the 2011 Tohoku earthquake.Earth,PlanetsandSpace, 66: 106.

Jacobs A, Sandwell D, Fialko Y. 2002. The 1999 (MW7.1) Hector Mine, California, earthquake: Near-field postseismic deformation from ERS interferometry.Bull.Seismol.Soc.Am., 92(4): 1433-1442. Jiang G Y, Xu C J, Wen Y M, et al. 2013. Inversion for coseismic slip distribution of the 2010MW6.9 Yushu Earthquake from InSAR data using angular dislocations.GeophysicalJournalInternational, 194(2): 1011-1022, doi: 10.1093/gji/ggt141. Jónsson S, Segall P, Pedersen R, et al. 2003. Post-earthquake ground movements correlated to pore-pressure transients.Nature, 424(6945): 179-183. Kanamori H. 1977. The energy release in great earthquakes.J.Geophys.Res., 82(20): 2981-2987.

Li C Y, Pang J Z, Zhang Z Q. 2012. Characteristics, geometry, and segmentation of the surface rupture associated with the 14 April 2010 Yushu Earthquake, Eastern Tibet, China.BulletinoftheSeismologicalSocietyofAmerica, 102(4): 1618-1638, doi: 10.1785/0120110261. Li Q. 2015. Coseismic deformation of the 2010MS7.1 Yushu earthquake [Master thesis] (in Chinese). Beijing: Institute of Earthquake Science, China Earthquake Administration, 47-51. Li Z H, Elliott J R, Feng W P, et al. 2011. The 2010MW6.8 Yushu (Qinghai, China) earthquake: Constraints provided by InSAR and body wave seismology.JournalofGeophysicalResearch, 116(B10): B10302, doi: 10.1029/2011JB008358.

Lin A M, Rao G, Jia D, et al. 2011. Co-seismic strike-slip surface rupture and displacement produced by the 2010Mw6.9 Yushu earthquake, China, and implications for Tibetan tectonics.JournalofGeodynamics, 52(3-4): 249-259.

Marone C J, Scholtz C H, Bilham R. 1991. On the mechanics of earthquake afterslip.J.Geophys.Res., 96(B5): 8441-8452. Massonnet D, Thatcher W, Vadon H. 1996. Detection of postseismic fault-zone collapse following the Landers earthquake.Nature, 382(6592): 612-616. Meng G J, Ren J W, Su X N, et al. 2013. Coseismic deformation of the 2010MW6.9 Yushu earthquake derived from GPS data.SeismologicalResearchLetters, 84(1): 57-64.

Nykolaishen L, Dragert H, Wang K L, et al. 2015. GPS observations of crustal deformation associated with the 2012MW7.8 Haida Gwaii earthquake.Bull.Seismol.Soc.Am., 105(2B): 1241-1252, doi: 10.1785/0120140177.

Okada Y. 1985. Surface deformation due to shear and tensile faults in a half-space.Bull.Seismol.Soc.Am., 75(4): 1135-1154.

Ozawa, S, Nishimura T, Suito H, et al. 2011. Coseismic and postseismic slip of the 2011 magnitude-9 Tohoku-Oki earthquake.Nature, 475(7356): 373-376, doi: 10.1038/nature10227. Peltzer G, Rosen P, Rogez F, et al. 1996. Postseismic rebound in fault step-overs caused by pore fluid flow.Science, 273(5279): 1202-1204.

Peltzer G, Rosen P, Rogez F, et al. 1998. Poroelastic rebound along the Landers 1992 earthquake surface rupture.J.Geophys.Res., 103(B12): 30131-30145. Perfettini H, Avouac J P. 2004. Postseismic relaxation driven by brittle creep: A possible mechanism to reconcile geodetic measurements and the decay rate of aftershocks, application to the Chi-Chi earthquake, Taiwan.J.Geophys.Res., 109(B2): B02304, doi: 10.1029/2003JB002488.

Pollitz F F. 1997. Gravitational viscoelastic postseismic relaxation on a layered spherical earth.J.Geophys.Res., 102(B8): 17921-17941. Qu C Y, Zhang G F, Shan X J, et al. 2013. Coseismic and postseismic deformation fields of the 2010 Yushu, QinghaiMS7.1 earthquake and their evolution processes.ChineseJ.Geophys. (in Chinese), 56(7): 2280-2291, doi: 10.6038/cjg20130715. Ren J W, Wang M. 2005. GPS measured crustal deformation of theMS8.1 Kunlun earthquake on November 14th 2001 in Qinghai-Xizhang Plateau.QuaternarySciences(in Chinese), 25(1): 34-44.

Rice J R, Cleary M P. 1976. Some basic stress diffusion solutions for fluid-saturated elastic porous media with compressible constituents.Rev.Geophys., 14(2): 227-241.

Ryder I, Wang H, Bie L D, et al. 2014. Geodetic imaging of late postseismic lower crustal flow in Tibet.EarthandPlanetaryScienceLetters, 404: 136-143.

Savage J C, Burford R O. 1973. Geodetic determination of relative plate motion in central California.JournalofGeophysicalResearch, 78(5): 832-845.

Savage J C, Svarc J L. 2009. Postseismic relaxation following the 1992M7.3 Landers and 1999M7.1 Hector Mine earthquakes, southern California.J.Geophys.Res., 114(B1): B01401, doi: 10.1029/2008JB005938. Segall P, Harris R. 1987. Earthquake deformation cycle on the San Andreas fault near Parkfield, California.J.Geophys.Res., 92(B10): 10511-10525.

Segall P. 2010. Earthquake and Volcano Deformation. Princeton: Princeton University Press, 321-324.

Shen Z K, Jackson D D, Feng Y J, et al. 1994. Postseismic deformation following the Landers earthquake, California, 28 June 1992.Bull.Seismol.Soc.Am., 84(3): 780-791.

Tobita M, Nishimura T, Kobaysahi T, et al. 2011. Estimation of coseismic deformation and a fault model of the 2010 Yushu earthquake using PALSAR interferometry data.EarthPlanet.Sci.Lett., 307(3-4): 430-438.

Wang K L, Hu Y, He J H. 2012. Deformation cycles of subduction earthquakes in a viscoelastic Earth.Nature, 484(7394): 327-332, doi: 10.1038/nature11032.

Wang L F, Hainzl S, Sinan ?zeren M, et al. 2010. Postseismic deformation induced by brittle rock damage of aftershocks.J.Geophys.Res., 115(B10): B10422, doi: 10.1029/2010JB007532. Wang R J, Lorenzo-Martin F, Roth F. 2006. PSGRN/PSCMP—a new code for calculating co-and post-seismic deformation, geoid and gravity changes based on the viscoelastic-gravitational dislocation theory.Computers&Geosciences, 32(4): 527-541.

Wang W L, Wu J P, Fang L H, et al. 2013. Relocation of the YushuMS7.1 earthquake and its aftershocks in 2010 from HypoDD.ScienceChinaEarthSciences, 56(2): 182-191, doi: 10.1007/s11430-012-4450-z.

Wang Y Z, Wang M, Shen Z K, et al. 2013. Inter-seismic deformation field of the Ganzi-Yushu fault before the 2010MW6.9 Yushu earthquake.Tectonophysics, 584: 138-143.

Wells D L, Coppersmith K J. 1994. New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement.Bull.Seismol.Soc.Am., 84(4): 974-1002. Wen X Z, Xu X W, Zheng R Z, et al. 2003. Average slip rate and present-day earthquake ruptures of the Ganzi-Yushu fault.

ScienceinChina,SeriesD(in Chinese), 33(supp.): 199-207.

Wen Y M, Li Z H, Xu C J, et al. 2012. Postseismic motion after the 2001MW7.8 Kokoxili earthquake in Tibet observed by InSAR time series.J.Geophys.Res., 117(B8): B08405, doi: 10.1029/2011JB009043.

Zhang Y, Xu L S, Chen Y T. 2010. Source process of the 2010 Yushu, Qinghai, earthquake.ScienceChinaEarthSciences, 53(9): 1249-1251.

附中文參考文獻(xiàn)

丁開華, 許才軍, 溫?fù)P茂. 2013. 汶川地震震后形變的GPS反演. 武漢大學(xué)學(xué)報·信息科學(xué)版, 38(2): 131-135.

李奇. 2015. 基于InSAR和GPS數(shù)據(jù)研究2010玉樹MS7.1級地震同震變形特征[碩士論文]. 北京: 中國地震局地震預(yù)測研究所, 47-51.

屈春燕, 張桂芳, 單新建等. 2013. 2010年青海玉樹地震同震-震后形變場特征及演化過程. 地球物理學(xué)報, 56(7): 2280-2291, doi: 10.6038/cjg20130715.

任金衛(wèi), 王敏. 2005. GPS觀測的2001年昆侖山口西MS8.1級地震地殼變形. 第四紀(jì)研究, 25(1): 34-44.

王未來, 吳建平, 房立華等. 2012. 2010年玉樹MS7.1地震及其余震的雙差定位研究. 中國科學(xué): 地球科學(xué), 42(7): 1037-1046.

聞學(xué)澤, 徐錫偉, 鄭榮章等. 2003. 甘孜—玉樹斷裂的平均滑動速率與近代大地震破裂. 中國科學(xué)(D輯), 33(增刊): 199-207.

張勇, 許力生, 陳運泰. 2010. 2010年青海玉樹地震震源過程. 中國科學(xué): 地球科學(xué), 40(7): 819-821.

(本文編輯 汪海英)

Postseismic deformation associated with the 2010 Yushu, Qinghai MS7.1 earthquake by GPS observations

MENG Guo-Jie1, SU Xiao-Ning1, XU Wan-Zhen1, REN Jin-Wei1, YANG Yong-Lin2, Nikolay V. Shestakov3,4

1KeyLaboratoryofEarthquakeForecast(InstituteofEarthquakeScience),ChinaEarthquakeAdministration,Beijing100036,China2InstituteofEngineeringSurveyandMapping,EarthquakeAdministrationofSichuanProvince,SichuanYa′an625000,China3InstituteofAppliedMathematics,FarEasternBranchofRussianAcademyofSciences,Vladivostok690091,Russia4FarEasternFederalUniversity,Vladivostok690068,Russia

The YushuMS7.1 earthquake occurred on April 14th, 2010 in the southeastern Tibetan Plateau, China. The seismic event ruptured the Yushu segment of the pre-existing Garzê-Yushu fault. We established 15 GPS stations over 7～10 days following the earthquake, which constituted a ～360 km long transection across the rupture trace, encompassing 1 continuously operating, 3 semi-continuous and 11 survey-mode sites. Most of the sites were observed for more than 240 days to capture the spatiotemporal changes of postseismic deformation. The secular deformation of the 15 sites was corrected using Euler vectors and a dislocation model. A layered visco-elastic dislocation model was employed to calculate aftershock-induced deformation, showing that aftershocks contributed considerably to the deformation of near-field sites. A logarithmic function was used to fit the post-seismic displacements, yielding a characteristic decay parameter of 6.7±1.2 days. The Steepest Decrease Method of optimization was employed to invert for spatiotemporal changes of afterslip. Afterslip was found mostly at the northwestern and southeastern flanks of coseismic slip. The afterslip patch at the northwestern flank almost reaches the surface, and the other one being directly beneath the major cosesimic slip patch. Maximum afterslip is located at 10～20 km depth range southeast of Gyêgu town. With the time lapse after the mainshock, the 2 domains of major afterslip were stationary in space, but their area was significantly increased. The moment released by the afterslip is (1.5～5.1)×1018Nm, equivalent to an earthquake ofMW6.1～6.4. Deformation due to poroelastic rebound cannot explain the observed surface deformation. The contribution to observed postseismic deformation from relaxation of lower crust and upper mantle is not significant, calculated using a Maxwell visco-elastic model. The rapid decay of postseismic deformation and the good fit of the after slip model suggest that the dominant source in the early stage of postseismic process is most likely ascribed to afterslip.

YushuMS7.1 earthquake; Postseismic deformation; Tectonic movement; Aftershock activity; Postseismic afterslip

10.6038/cjg20161219.

國家自然科學(xué)基金(4117004，41461164004)，國家國際科技合作項目(2015DFR21100，2010DFB20190)，2010年玉樹7.1級地震科學(xué)考察(0210660),中國地震局地震預(yù)測研究所基本科研業(yè)務(wù)費專項重點項目(2014IES010102)，俄羅斯遠(yuǎn)東聯(lián)邦大學(xué)項目(14-08-01-05_m)聯(lián)合資助.

孟國杰，男，1969年生，研究員，主要研究方向為地震形變，構(gòu)造活動，地球動力學(xué).E-mail:guojiemeng@126.com

10.6038/cjg20161219

P315

2015-10-15，2016-10-28收修定稿

孟國杰，蘇小寧，徐婉楨等. 2016. 基于GPS觀測研究2010年青海玉樹MS7.1地震震后地殼形變特征及其機(jī)制. 地球物理學(xué)報，59(12):4570-4583,

Meng G J, Su X N, Xu W Z,et al. 2016. Postseismic deformation associated with the 2010 Yushu,QinghaiMS7.1 earthquake by GPS observations.ChineseJ.Geophys. (in Chinese),59(12):4570-4583,doi:10.6038/cjg20161219.