基于非均勻地殼模型的格林函數(shù)法反演2001年昆侖山口西Ms8.1地震的同震滑移*

周 宇 聶兆生 賈治革 陳正松 楊少敏

(中國(guó)地震局地震研究所(地震大地測(cè)量實(shí)驗(yàn)室)，武漢 430071)

基于非均勻地殼模型的格林函數(shù)法反演2001年昆侖山口西Ms8.1地震的同震滑移*

周 宇 聶兆生 賈治革 陳正松 楊少敏

(中國(guó)地震局地震研究所(地震大地測(cè)量實(shí)驗(yàn)室)，武漢 430071)

建立基于地殼不均勻的三維有限元模型，將昆侖山口西Ms8.1地震產(chǎn)生的長(zhǎng)約400 km的破裂帶在走向上分為16段，生成各個(gè)破裂段滑移相對(duì)于地表位移的格林函數(shù)，并以GPS觀(guān)測(cè)的水平位移為約束，用阻尼最小二乘法反演破裂滑移分布，結(jié)果表明反演結(jié)果與地質(zhì)考察的地表破裂數(shù)據(jù)接近。用反演結(jié)果重構(gòu)的地表位移場(chǎng)與觀(guān)測(cè)值一致，并體現(xiàn)了斷裂帶南北兩側(cè)不同的形變模式。

有限元;最小二乘法;昆侖山口西Ms8.1地震;同震滑移;GPS

1 引言

昆侖山口西Ms8.1地震發(fā)生在東昆侖斷裂帶上，該斷裂帶是青藏高原一條古老的縫合帶，其左旋走滑非?；钴S。震后，美國(guó)哈佛大學(xué)、美國(guó)地質(zhì)調(diào)查局(USGS)、日本東京大學(xué)地震研究所和中國(guó)地震局地球物理研究所先后給出了此次地震的震中位置及矩張量震源機(jī)制解，結(jié)果略有差異，但多認(rèn)為此次主震的破裂方式以走滑破裂為主，兼有少量的逆沖滑動(dòng)分量［1，2］。野外地質(zhì)考察發(fā)現(xiàn)，地震沿東昆侖斷裂帶西段形成了長(zhǎng)約350～420 km的地表破裂帶，位錯(cuò)以走滑為主，最大左旋水平錯(cuò)距 6～7 m［3-5］。

文獻(xiàn)［1，2］運(yùn)用全球地震波資料反演了該地震的時(shí)空破裂過(guò)程。文獻(xiàn)［1］的反演結(jié)果為:斷層的破裂長(zhǎng)度約380 km、寬約30 km、最大滑動(dòng)量5.8 m、平均錯(cuò)動(dòng)(左旋走滑)1.8 m;文獻(xiàn)［2］假定斷層破裂深度為40 km，得到的斷層最大滑動(dòng)量為2.2 m，平均滑動(dòng)量?jī)H1.2 m。顯然這些結(jié)果與地表破裂的實(shí)測(cè)滑動(dòng)量［3-5］差別較大。

為使反演結(jié)果更加真實(shí)地表現(xiàn)地表破裂現(xiàn)象及斷層的活動(dòng)狀況，文獻(xiàn)［6～10］分別以大地測(cè)量結(jié)果為約束，計(jì)算了該地震的同震位移，并討論了地表形變特征;單新建［11］則利用InSAR資料給出了衛(wèi)星視線(xiàn)向的同震形變場(chǎng)。無(wú)論是用GPS還是InSAR，地表形變的觀(guān)測(cè)結(jié)果都表明，東昆侖斷裂帶南北兩側(cè)的形變顯示出明顯的非對(duì)稱(chēng)性，斷層南側(cè)的同震位移和水平位移梯度的南北向分量都明顯大于斷層北側(cè)，表明斷層南北兩側(cè)介質(zhì)的彈性力學(xué)性質(zhì)有明顯的差異。這些反演結(jié)果與實(shí)地考察結(jié)果具有較好的一致性，地震波速資料也證實(shí)了反演結(jié)果的可信性［12-14］。

但這些研究使用的是Okada模型的解析公式，而Okada模型基于均勻彈性半空間，所以無(wú)法解釋南北盤(pán)形變特征的明顯差異。要解決這一問(wèn)題，必須加入介質(zhì)在橫向和垂向上的不均勻性，然而目前的位錯(cuò)理論并沒(méi)有非均勻介質(zhì)條件下的解析解，所以采用數(shù)值模擬方法就成為必然。王輝等［15］用三維有限元數(shù)值方法，考慮介質(zhì)的不均勻性模擬了昆侖山口西地震的同震形變場(chǎng)，但該研究的目的在于重構(gòu)空間連續(xù)的形變場(chǎng)，將野外地質(zhì)勘查獲得的地表走滑量作為地震破裂分布。然而，并非所有的破裂都會(huì)到達(dá)地表，而地表的破裂不一定和深部的破裂一致。Masterlark［16］發(fā)現(xiàn)，在板塊俯沖帶同震及震后形變的數(shù)值模擬中，在非均勻地殼模型中用基于Okada模型反演得到的滑動(dòng)分布來(lái)模擬地表形變場(chǎng)，可能會(huì)導(dǎo)致不可忽視的誤差。為解決這個(gè)問(wèn)題，本文提出用三維有限元建立非均勻的介質(zhì)模型，計(jì)算各個(gè)破裂段對(duì)于地表位移的格林函數(shù)，然后用GPS觀(guān)測(cè)數(shù)據(jù)反演各個(gè)破裂段的滑移。

2 模型的建立

昆侖山地震后，中國(guó)地震局地震研究所聯(lián)合中國(guó)地震局第二地形變監(jiān)測(cè)中心迅速開(kāi)展了應(yīng)急觀(guān)測(cè)，獲得了40個(gè)點(diǎn)的同震位移數(shù)據(jù)。我們將利用這些數(shù)據(jù)，以單日時(shí)段為基本單位解算各GPS點(diǎn)的三維坐標(biāo)。由于地震位錯(cuò)以走滑為主，GPS站點(diǎn)的垂直位移不大，而GPS觀(guān)測(cè)在垂向上的誤差比水平分量的要高出一個(gè)數(shù)量級(jí)，所以本文拋棄垂向位移，只用水平位移。

2.1 單元?jiǎng)澐?/h3>

考慮巖石圈縱向上的分層和斷裂帶南北兩側(cè)的分界，我們依據(jù)亞?wèn)|格爾木額濟(jì)納旗地學(xué)斷面建立的巖石圈地質(zhì)模型［9］，并參考此斷面得到的介質(zhì)密度、速度、泊松比和莫霍面深度，用縱波速和泊松比換算出楊氏彈性模量。研究區(qū)介質(zhì)的分區(qū)見(jiàn)圖1，參數(shù)的設(shè)置見(jiàn)表1，單元的劃分如圖2。

圖1 模型剖面示意圖Fig.1 Sketch of model profile

表1 介質(zhì)參數(shù)Tab.1 Medium parameters

模型中我們將斷裂同一側(cè)的介質(zhì)視為橫向均勻。

2.2 三維有限元網(wǎng)格與邊界條件

網(wǎng)格總體長(zhǎng)1 500 km，寬1 000 km，深200 km，頂面的中心與斷裂帶中點(diǎn)重合。斷裂帶用長(zhǎng)400 km、深14 km，完全與地面垂直的矩形面表示，走向277°(圖2)。模型的網(wǎng)格采用四面體單元，破裂面附近形變梯度大的地方做精細(xì)劃分，遠(yuǎn)場(chǎng)位移梯度小的地方做粗疏劃分。數(shù)據(jù)試驗(yàn)表明，在破裂面附近的網(wǎng)格單元的邊長(zhǎng)應(yīng)該比破裂面的邊長(zhǎng)小一個(gè)數(shù)量級(jí)，才能保證有限元數(shù)值計(jì)算的最大誤差在2%以?xún)?nèi)，因此在斷裂帶附近的網(wǎng)格邊長(zhǎng)設(shè)定為2～4 km，在邊界處邊長(zhǎng)為35～50 km。網(wǎng)格共有節(jié)點(diǎn)58 349個(gè)，四面體單元313 502個(gè)(圖2)。

約束條件為:網(wǎng)格4個(gè)側(cè)面的邊界的水平位移為零，底面的垂向位移為零，側(cè)面和底面其他方向上的位移分量自由，地表位移的3個(gè)分量均為自由。

使用Pylith軟件進(jìn)行有限元計(jì)算。

圖2 三維有限元網(wǎng)格的水平面Fig.2 Horizontal surface of 3D-FEM mesh

2.3 格林函數(shù)的生成

設(shè)任意地面點(diǎn)xi的位移矢量vi為關(guān)于斷層某個(gè)區(qū)域ξj滑移矢量uj的線(xiàn)性函數(shù):

式中g(shù)(u)為格林函數(shù)。

將式(1)寫(xiě)成矩陣形式，有:

由于觀(guān)測(cè)條件所限，GPS點(diǎn)不多，分布也不夠理想，因此沒(méi)有基于破裂模式在平面上對(duì)斷裂帶做不等長(zhǎng)、垂向上做不等深的精細(xì)劃分，為方便起見(jiàn)，將破裂帶分為16個(gè)等長(zhǎng)的分段，每段為長(zhǎng)25 km、深14 km的矩形。并設(shè)定滑動(dòng)向量中的一個(gè)分量為1，其他分量為零。用三維有限元計(jì)算正演各GPS站的位移。將2M個(gè)位移向量從左到右并列即可獲得模型的格林函數(shù)矩陣。

3 反演方法

以最小方差為反演目標(biāo)函數(shù)，用廣義逆反演超定方程(2)，有:

其中CV是觀(guān)測(cè)向量V的協(xié)方差矩陣。

由于各GPS站觀(guān)測(cè)值獨(dú)立，且同一個(gè)GPS站的位移觀(guān)測(cè)值的兩個(gè)水平分量相關(guān)系數(shù)接近于零，故CV為對(duì)角陣，有:

σi是觀(guān)測(cè)向量第i個(gè)分量的標(biāo)準(zhǔn)差。為表示方便，令:

代入方程(3)有

模擬結(jié)果表明，用方程(8)反演，其結(jié)果發(fā)生了畸變，破裂面滑移的絕對(duì)值遠(yuǎn)遠(yuǎn)超出了正常范圍。經(jīng)分析，發(fā)生畸變是因?yàn)榫仃嘒'TG'的最大和最小特征值相差太大，比率高達(dá)1011，導(dǎo)致反演的病態(tài)，所以無(wú)法用解決超定問(wèn)題的方法來(lái)反演。由此本文改用阻尼最小二乘法進(jìn)行反演，由于阻尼最小二乘法加入阻尼因子避免了病態(tài)問(wèn)題。反演方程如下:

其中:β＞0，為阻尼因子;I為單位矩陣。殘差大小為:

其中norm()為向量的2-范數(shù)，tr()為矩陣的跡。

反演結(jié)果顯示，雖然?U相差很大，但都在GPS點(diǎn)上產(chǎn)生了幾乎一樣的位移。分析其原因，一是觀(guān)測(cè)數(shù)據(jù)不足，導(dǎo)致數(shù)據(jù)向量對(duì)模型向量不敏感;二是阻尼因子太小無(wú)法解決反演的病態(tài)問(wèn)題。由此我們引入地震矩來(lái)選取合適的阻尼因子。地震矩公式如下:

其中Aj是第j個(gè)破裂面的面積，sj是第j個(gè)破裂面的滑動(dòng)距離，μ是剪切模量。哈佛矩心張量的地震矩為5.9×1020Nm，與之對(duì)應(yīng)的阻尼因子β=10-3.5。

最終，反演獲得的破裂面滑動(dòng)分布的最優(yōu)解如圖3所示。

圖3 基于非均勻地殼模型的同震滑移分布的反演結(jié)果Fig.3 Inversed coseismic slip distribution based on the heterogeneous crust model

4 反演結(jié)果分析

從圖3可以看到，反演得到的同震滑移以走滑為主，傾滑分量在大多數(shù)分段上小于50 cm。左旋走滑分布有三個(gè)峰值，從大到小分別分布于東經(jīng)93.6°、92.6°與94.4°處，其量級(jí)分別為5.45、4.84和3.82 m，與地質(zhì)考察的結(jié)果［5］基本一致。但傾滑分量在兩個(gè)分段上超過(guò)了1 m，與實(shí)際不符。

將反演獲得的滑動(dòng)分布，代入公式(2)計(jì)算得到各GPS點(diǎn)的水平位移矢量如圖4所示。

圖4 水平位移的模擬值與GPS觀(guān)測(cè)值的比較Fig.4 Comparison between the simulated horizontal displacement field and the GPS observations

圖4中，紅色箭頭是GPS觀(guān)測(cè)得到的水平位移，綠色箭頭是模擬的水平位移，誤差橢圓置信區(qū)間為95%。其中BS33和WT02兩個(gè)點(diǎn)的位移超出圖幅，為表示方便，圖中比例尺所示的位移大小為真實(shí)值的一半。在大部分觀(guān)測(cè)點(diǎn)上，觀(guān)測(cè)值和模擬值都比較一致，在幾個(gè)離斷層較近的點(diǎn)上二者幾乎相同。模擬位移也很好地體現(xiàn)了斷裂帶南側(cè)的水平位移的南北向梯度比北側(cè)的大。

BS26、G0CQ兩個(gè)點(diǎn)的震前觀(guān)測(cè)時(shí)間分別在1995年和1997年［7］，而1995—2001年該地區(qū)發(fā)生5級(jí)以上地震10次，所以GPS觀(guān)測(cè)位移可能包含了這些地震引起的位移，使得觀(guān)測(cè)值遠(yuǎn)大于模擬值。

在斷裂南側(cè)，除了近場(chǎng)的WT02，其他幾個(gè)點(diǎn)的模擬位移比實(shí)測(cè)的位移小，但在各點(diǎn)上二者相差的比例幾乎相同，這說(shuō)明模擬結(jié)果的位移梯度或者應(yīng)變場(chǎng)與觀(guān)測(cè)結(jié)果幾乎是相同的。引起位移不一致的原因可能是在WT02和WUDA之間存在著彈性模量較大的局部區(qū)域，導(dǎo)致該區(qū)域的應(yīng)變較小，使得從WT02到WUDA的位移梯度變小。

處于阿爾金斷裂帶附近的幾個(gè)點(diǎn)上(圖4中的西北角)，模擬位移和觀(guān)測(cè)位移的大小比較接近，但是觀(guān)測(cè)位移為北東向，模擬位移為北西向，二者近于垂直。經(jīng)分析可能是沒(méi)有用軟弱帶模擬阿爾金斷裂所致。

從整體上看，斷裂帶西側(cè)的模擬值與觀(guān)測(cè)值的差別明顯大于東側(cè)，從反演的角度來(lái)看，可能是斷裂帶西半段的近場(chǎng)沒(méi)有觀(guān)測(cè)點(diǎn)約束之故。

圖5 基于彈性半空間模型的同震滑移分布的反演結(jié)果Fig.5 Inversed coseismic slip distribution based on the homogeneous elastic half space model

圖6 非均勻地殼模型對(duì)于兩種不同模型的滑動(dòng)反演結(jié)果的地表響應(yīng)Fig.6 Ground surface response of heterogeneous crust model to two kinds of inversed slip distribution

GPS觀(guān)測(cè)不能提供同震位移場(chǎng)在空間上的連續(xù)分布，在反演獲得滑移分布之后帶入有限元模型正演，可以重建全區(qū)域的位移場(chǎng)和應(yīng)變場(chǎng)。為了對(duì)比，本文建立了均勻彈性半空間的模型，用Okada公式計(jì)算格林函數(shù)，并用同樣的反演方法計(jì)算了同震滑移分布。雖然均勻彈性半空間模型也可以在觀(guān)測(cè)點(diǎn)上產(chǎn)生與非均勻地殼模型幾乎相同的模擬位移，然而其反演的滑移分布(圖5)，與非均勻地殼模型有明顯差別，雖然起伏趨勢(shì)基本一致，但其量值有差別，如走滑分量最大值為7.2 m，最小值為-1.5 m，傾滑分量最大值為3.6 m，并出現(xiàn)了右旋走滑，且傾滑分量過(guò)大，顯然不如非均勻模型更接近真實(shí)情形。將均勻彈性半空間模型的滑移反演結(jié)果代入非均勻模型中，兩種模型的位移場(chǎng)比較如圖6所示，在絕大多數(shù)點(diǎn)上，均勻彈性半空間模型的位移比非均勻模型的位移小30%以上，這說(shuō)明基于地殼非均勻性反演同震滑移更加合理，而在非均勻性模型中直接使用均勻模型反演得到滑移分布是值得商榷的。

6 結(jié)論

1)阻尼最小二乘法可以克服反演中的病態(tài)問(wèn)題，而不必事先限定滑移分布的上下界。地質(zhì)考察結(jié)果驗(yàn)證了反演結(jié)果的可靠性。雖然沒(méi)有水準(zhǔn)測(cè)量數(shù)據(jù)，使得斷裂傾滑分量的反演結(jié)果不夠準(zhǔn)確，然而在觀(guān)測(cè)數(shù)據(jù)如此少的情況下能得到這樣的結(jié)果，充分證明了反演方法的有效性。

2)基于地殼橫向和縱向不均勻的三維有限元模型模擬的同震位移場(chǎng)與觀(guān)測(cè)值基本一致，雖然由于模型沒(méi)有考慮地殼東西向的不均勻性和其他局部特征，導(dǎo)致在少數(shù)觀(guān)測(cè)點(diǎn)上與觀(guān)測(cè)值差別略大，但還是很好地解釋了東昆侖斷裂帶南北兩側(cè)形變模式的不同。

3)在沒(méi)有地表滑移分布的情況下，如果地殼存在明顯的橫向不均勻性，將均勻彈性半空間模型的滑移反演結(jié)果作為內(nèi)邊界條件代入非均勻地殼模型來(lái)模擬形變場(chǎng)是不恰當(dāng)?shù)摹?/p>

1 Lin A，Kikuchi M and Fu B.Rupture segmentation and process of the 2001 Mw7.8 central Kunlun earthquake，China［J］.Bull Seism Soc Amer.，2003，93:2 477-2 492.

2 許力生，陳運(yùn)泰.用長(zhǎng)周期波形資料反演2001年11月14日昆侖山口地震時(shí)空破裂過(guò)程［J］.中國(guó)科學(xué)(D輯)，2004，34(3):256-264.(Xu Lisheng and Chen Yuntai.Temporal and spatial rupture process of the great Kunlun Mountain Pass earthquake of Nov.14，2001，from GDSN long period waveform data［J］.Science in China(Ser D)，2004，34(3):256-264)

3 徐錫偉，等.2001年11月14日昆侖山庫(kù)湖地震(Ms8.1)地表破裂帶的基本特征［J］.地震地質(zhì)，2002，24(1):1-13.(Xu Xiwei，et al.Characteristic featur es of the surface ruptures of the HohSai Hu(Kunlunshan)earthquake(Ms 8.1)northern Tibetan Plateau，China［J］.Seismology and Geology，2002，24(1):1-13)

4 Lin Aiming，et al.Coseismic strikeslip and rupture length produced by the 2001 Ms8.1 central Kunlun earthquake［J］.Science，2002，296:2 015-2 017.

5 陳杰，等.2001年昆侖山口西Ms8.1地震地表同震位移分布特征［J］.地震地質(zhì)，2004，26(3):378-392.(Chen Jie，et al.Surface rupture zones of the 2001 earthquake Ms 8.1 west of Kunlun Pass，northern Qinghai-Tibet plateau［J］.Quaternary Sciences，2003，23(6):629-639)

6 喬學(xué)軍，等.昆侖山口西Ms8.1地震的地殼變形特征［J］.大地測(cè)量與地球動(dòng)力學(xué)，2002，(4):6-11.(Qiao Xuejun，et al.Characteristics of crustal deformation relating to Ms8.1 Kunlunshan earthquake［J］.Journal of Geodesy and Geodynamics，2002，(4):6-11)

7 任金衛(wèi)，王敏.GPS觀(guān)測(cè)的2001年昆侖山口西Ms8.1地震地殼變形［J］.第四紀(jì)研究，2005，25(1):34-44.(Ren Jinwei and Wang Min.GPS measured crustal deformation of the Ms8.1 Kunlun earthquake on November 14th 2001 in Qinhai-Tibet Plateau［J］.Quaternary Sciences，2005，25 (1):34-44)

8 譚凱，王琪，申重陽(yáng).用大地測(cè)量數(shù)據(jù)反演2001年昆侖山地震［J］.大地測(cè)量與地球動(dòng)力學(xué)，2004，(3):47-50.(Tan Kai，Wang Qi and Shen Chongyang.Using geodetic data to inverse coseismic dislocation of 2001 Kunlun earthquake［J］.Journal of Geodesy and Geodynamics，2004，(3): 47-50)

9 萬(wàn)永革，等.利用GPS和水準(zhǔn)測(cè)量資料反演2001年昆侖山口西8.1級(jí)地震的同震滑動(dòng)分布［J］.地震地質(zhì)，2004，26(3):393-404.(Wan Yongge，et al.Co-seismic slip distribution of the 2001 west of kunlun mountain pass earthquake inverted by GPS and levelling data［J］.Seismology and Geology，2004，26(3):393-404)

10 馬超，單新建.昆侖山Ms8.1地震震源參數(shù)的多破裂段模擬研究［J］.地球物理學(xué)報(bào)，2006，49(2):428-437.(Ma Chao and Shan Xinjian.A multi-segment analytic modeling of hypocentral geometric characteristic parameters of the Ms8.1 earthquake at the Kunlun mountains［J］.Chinese J Geophys.，2006，49(2):428-437)

11 單新建，柳稼航，馬超.2001年昆侖山口西8.1級(jí)地震同震形變場(chǎng)特征的初步分析［J］.地震學(xué)報(bào)，2004，26 (5):474-480.(Shan Xinjian，Liu Jiahang and Ma Chao.Preliminary analysis on characteristics of coseismic deformation associated with Ms8.1 western Kunlunshan Pass earthquake in 2001［J］.Acta Seimologica Sinica，2004，26 (5):474-480)

12 吳功建，高銳，欽范.青藏高原“亞?wèn)|-格爾木地學(xué)斷面”綜合地球物理調(diào)查與研究［J］.地球物理學(xué)報(bào)，1991，34 (5):552-561.(Wu Gongjian，Gao Rui and Yu Qinfan.Integrated inverstigation of the Qinghai-Tibet plateau along the Yadong-Golmud geosciences transect［J］.Chinese J Geophys.，1991，34(5):552-562)

13 Zhu L and Helmberger D V.Moho offset across the northern margin of the Tibetan plateau［J］.Science，1998，281: 1 170-1 172.

14 Okada Y.Surface deformation due to shear and tensile faults in a half space［J］.Bull seism Soc Am.，1985，75: 1 135-1 154.

15 王輝，等.昆侖山口西Ms8.1地震同震影響場(chǎng)的數(shù)值模擬［J］.地震地質(zhì)，2007，29(3):637-647.(Wang Hui，et al.The numerical simulation on coseismic effect of the 14 November 2001 great Kunlun earthquake，Northern Tibet，China［J］.Seismology and Geology，2007，29(3):637-647)

16 Masterlark T，Demets C and Wang H F.Homogeneousvs heterogeneous subduction zone models:Coseismic and postseismic deformation［J］.J Geophys Res.，2001，28(21): 4 047-4 050.

INVERSION OF COSEISMIC SLIP DISTRIBUTION OF 2001 Ms8.1 KUNLUN EARTHQUAKE FROM GPS OBSERVATIONS WITH GREEN’S FUNCTION METHOD BASED ON A HETEROGENEOUS CRUST MODEL

Zhou Yu，Nie Zhaosheng，Jia Zhige，Chen Zhengsong and Yang Shaomin
(Key Laboratory of Earthquake Geodesy，Institute of Seismology，CEA，Wuhan 430071)

2001 Kunlun Ms8.1 earthquake ruptured as long as about 400 km.The displacement gradient perpendicular to the strike of rupture belt is significantly greater on the south side than that on the north.This fact cannot be explained by Okada model which is based on elastically homogeneous half space.A horizontally and vertically heterogeneous 3D-FEM model is established to generate Green’s functions of slip on 16 segments of the rupture belt with respect to the ground displacement field.With the damping LS(least square)algorithm，the slip distribution of rupture is inversed from observed horizontal displacements of 40 GPS sites，and the result is in good agreement with that from geological field investigation.The ground displacement field rebuilt by the inversed result is close to that observed by GPS，demonstrating different deformation patterns on the two sides of the rupture.

finite element;least square algorithm;Kunlun Ms8.1 earthquake;coseismic slip;GPS

1671-5942(2012)03-0022-05

2012-02-24

中國(guó)地震局地震行業(yè)科研專(zhuān)項(xiàng)(201208009)

周宇，男，1983年生，碩士，研究實(shí)習(xí)員，主要從事地震形變與活動(dòng)構(gòu)造的數(shù)值模擬與動(dòng)力學(xué)研究.E-mail:zhouyuking@foxmail.com

P315.72+5

A