GPS約束下唐山地區斷裂帶現今運動速度特征非連續接觸模擬*

孟慶筱

1)中國地震局第一監測中心，天津 300180

2)中國地震局地震研究所(地震大地測量重點實驗室)，武漢 430071

3)中國地震局地殼應力研究所武漢科技創新基地，武漢 430071

GPS約束下唐山地區斷裂帶現今運動速度特征非連續接觸模擬*

孟慶筱1，2，3)

1)中國地震局第一監測中心，天津 300180

2)中國地震局地震研究所(地震大地測量重點實驗室)，武漢 430071

3)中國地震局地殼應力研究所武漢科技創新基地，武漢 430071

以1999～2007年GPS觀測結果為基礎，使用二維接觸有限元方法計算唐山地區主要斷裂帶的運動速度狀態。模擬結果表明，現今地殼形變運動狀態下唐山地區主要斷裂帶運動狀態具有差異性。由于山海關塊隆對唐山塊陷的推擠作用，唐山斷裂帶以右旋拉張運動為主，同時在薊運河斷裂帶位置處保持著快速的擠壓應力積累，具有較高的地震危險性。

唐山;GPS;有限元;接觸單元;斷裂帶

研究認為，1976年唐山7.8級強烈地震與深大斷裂沒有表現出直接聯系，而是發生在唐山塊陷內部北東向的唐山斷裂帶上。邊界斷裂的存在使得被圍陷的唐山塊體得到暫時的平衡，形成相對鎖閉的區域。區域構造應力場發生急劇變化時，鎖閉平衡狀態被打破，從而發生唐山大地震［1］。2012-05-20唐山地區又發生M4.7地震，說明唐山斷裂構造的活動性明顯增強。研究唐山地區斷裂帶的現今形變特征，對于唐山地區地震危險性評價具有重要意義。本文通過構建唐山地區二維有限元接觸模型，利用中國大陸網絡工程復測資料1999～2007年GPS結果進行邊界條件約束，使用經過信度檢驗的有限元模擬結果分析唐山斷裂、豐臺-野雞坨斷裂、薊運河斷裂、寧河-昌黎斷裂、大李莊-古馬斷裂等唐山地區主要斷裂帶的現今形變特征，為唐山地區地震危險性分析提供依據。

1 唐山地區構造概況

唐山塊陷位于陰山-燕山斷褶帶和冀魯斷塊坳陷之間，四周被深大斷裂所包圍，塊體中央有唐山斷裂帶，唐山地震就發生在這個塊體內部的斷裂帶上。唐山塊陷在形態上為一個走向北東東的菱形塊體，四周被具有新活動的4條大斷裂和深斷裂所包圍，斷塊的南界為寧河-昌黎深斷裂，北界為豐臺-野雞坨大斷裂，二者走向均為北東東;東界是灤縣-樂亭大斷裂，西界為薊運河深斷裂，二者走向均成北西向。其活動特征分述如下［1］:

1)寧河-昌黎深斷裂:為唐山塊陷和樂亭塊陷的分界線，向北東延入渤海，向西與滄東斷裂相接。斷裂產生于呂梁運動，新活動強烈，斷裂兩盤第四紀沉積物落差500 m以上。在斷裂帶的西北盤上，有一平行于該斷裂的北東東向河流縱比降轉折帶，與斷裂的北西盤上升、南東盤下降有關。

2)豐臺-野雞坨大斷裂:為薊縣塊隆和唐山塊陷的交界線，控制古生代地層的發育。斷裂以北為巨厚的震旦系地層，以南主要為古生代地層，西段斷面傾向西北，東段傾向南東。在榛子鎮附近為斷面轉折的樞紐部位，新構造運動明顯。東段是山區和平原的分界線，西段斷裂兩側第四系厚度之差可達300～400 m，豐潤西南沿斷裂出現第四系沉積物厚度等值線的密集帶。

3)灤縣-樂亭大斷裂:斷裂控制山海關塊隆和唐山塊陷的活動，長期前者上升、后者下沉。斷面傾向北東，為高角度的逆沖斷層。在馬城-野雞坨之間有4條平行斷層，形成于前古生代，新構造運動強烈。

4)薊運河深斷裂:為唐山塊陷和滄縣塊隆的分界線，向西北延伸，與寶坻-香河東西向斷裂相連。香河寶坻段為重力梯度帶，東南段梯度緩斷裂控制古生代地層的發育和中生代盆地的形成。斷裂形成于前古生代，新構造運動強烈，不但控制薊運河的流向，而且以斷裂為界，兩側水系流向截然不同。

5)唐山斷裂帶:該斷裂由唐山-古冶斷裂、陡河斷裂和唐山-巍山-長山南坡斷裂3條平行斷裂組成，其中唐山-古冶斷裂南西段走向N30°E，東北段走向N50°E;陡河斷裂東北段為一向北西傾的正斷層，西南段由平行的4條小斷層組成，斷面傾向北西。最西一條為正斷層，東邊3條為逆斷層，全長約50 km;唐山-巍山-長山南坡斷裂由一些斷斷續續的北東向斷層組成，多為向北西傾的逆斷層。為便于建模，將以上3條斷裂合并為唐山斷裂帶。

根據以上論述，多期地質構造運動導致唐山地區呈現為被深大斷裂圍陷、分割的復雜塊體結構。由于地殼深部各圈層之間存在的相互作用在構造線位置最為突出，因此斷裂帶的展布對于地殼運動的水平速度場具有明確的控制作用，需要建立足夠真實的有限元模型來反映研究區內斷裂構造的復雜性以及相鄰塊體之間沿斷裂可能發生的相對滑動。

2 模型和方法

2.1 唐山地區二維有限元接觸模型

研究范圍為 117.5 ～ 119.0°E、38.9 ～ 40.4°N內的矩形區域。根據鄧啟東［2］給出的唐山地區構造情況和斷裂帶展布，本文構建了1個二維矩形彈性薄片模型(圖2)，主要涵蓋Ⅰ-冀中塊陷、Ⅱ-樂亭塊陷、Ⅲ-山海關塊隆、Ⅳ-薊縣塊隆、Ⅴ-唐山塊陷(圖1(a))。本研究使用GPS觀測結果作為邊界條件進行約束，其時間區段為1999～2007年。由于模擬時間尺度較短且運動速率較小，地殼變形過程基本處于穩定狀態。因此，在本研究中忽略動力學效應，將模型簡化為靜力學問題。本文采用二維彈性有限元模型模擬地殼水平運動，可以大大降低計算成本，將其作為精確三維地質體模型的合理近似，從而實現淺層巖石圈模擬［3］。

為了提高模型的可分析性，在保證模型客觀反映基本問題的前提下，對研究區內的次級斷裂帶進行簡化。納入模型的斷層共8條:(1)豐臺-野雞坨斷裂、(2)薊運河斷裂、(3)唐山斷裂、(4)寧河-昌黎斷裂、(5)灤縣-樂亭斷裂、(6)滄西斷裂、(7)灤縣-盧龍斷裂、(8)大李莊-古馬斷裂(圖1(a))。因此，該模型可以滿足唐山地區斷裂帶現今變形模擬分析的需求。

采用ANSYS有限元模擬環境，模型中構造塊體和易變形帶采用彈性連續本構關系，采用面單元對研究區進行單元剖分，見表1［4］。斷裂帶作為本研究中的重要因素，為了更加真實地反映斷裂帶附近破碎巖體較圍巖強度更低的實際情況，本文通過介質物性參數弱化的方法在模型中構建了寬度約為斷層長度1%的斷裂過程帶(圖3)［5］。

圖1 唐山地區地震構造、GPS站速度和有限元模擬結果Fig.1 Seismic tectonics，velocity by GPS observation and FEM simulation results in Tangshan area

圖2 有限元模型Fig.2 FEM model

圖3 斷裂過程帶示意圖［5］Fig.3 Schech of the process zone

為了充分考慮相鄰塊體可能沿著斷裂發生的相對滑動，將斷裂帶作為不連續面進行處理。使用conta175、targe169并以接觸對的形式進行建模，即接觸面之間用假想的法向和切向彈簧鏈接，其力學性質通過法向接觸剛度、切向接觸剛度和滑動摩擦系數來控制(表 2)［2，6］。

表1 Plane182單元的物理力學參數［4］Tab.1 Mechanical parameters of the unit no.182

表2 非接觸單元的物理力學參數［2］Tab.2 Mechanical parameters of non-contact unit

2.2 GPS數據及邊界條件

采用ITRF框架下中國大陸1999～2007年站速度結果，通過坐標旋轉，得到區域無旋轉坐標框架下GPS速度場(圖1(b))。該結果反映了唐山地區地殼的水平運動狀態，并可以清晰地觀察斷層兩側以及不同塊體上的GPS速度特征。將靠近模型矩形區域邊界處的GPS測點的實際觀測速度作為距離最近處節點的位移約束，通過樣條插值方法得到關鍵控制點之間其他節點的邊界條件。

3 計算結果與討論

3.1 地殼運動模擬結果及其信度檢驗

使用ANSYS計算得到的地殼運動模擬速度與GPS實際觀測結果的對比如圖1(b)。可以發現，GPS觀測值和有限元模擬值的一致性較好，說明模擬結果可以較好地反映研究區地殼運動的實際。圖4(a)、(b)分別為模擬結果以及與GPS觀測結果的東西向分量、南北向分量的誤差比對結果。可見，大部分測站的模擬速度都在GPS觀測誤差范圍之內。東西方向分量和南北方向分量的總體標準偏差分別為1.52和1.34 mm/a，模擬結果信度較高，可以直接用于唐山地區斷裂帶現今形變特征研究。

圖1(b)的地殼運動模擬結果顯示，唐山地區地殼形變場具有較為明顯的不均勻性和分界性，總體呈現出向西運動的速度場狀態，研究區自北向南，除個別點外，速度矢量方向由正西向逐漸轉變為南西向，不均勻性的同時，呈現較為規律的扇形分布。唐山塊陷內部GPS觀測點運動速率明顯小于樂亭塊陷、山海關塊隆和薊縣塊隆等外圍地塊，說明唐山塊陷在四周深大斷裂的圍陷下，在1999～2007年間呈現為相對鎖閉平衡狀態。

3.2 斷裂帶運動速度分析

依據單元接觸對的位移可以得到斷層運動速率，包括平行于斷層的水平走滑分量(右旋/左旋，右旋為正)和垂直斷層面的水平拉張分量(拉張/擠壓，拉張為正)，見表3和圖5。

本文選擇5條斷裂帶進行重點分析。運動速率計算結果顯示:

1)豐臺-野雞坨斷裂:該斷裂主要以左旋擠壓運動為主，其中左旋走滑速率為2.61 mm/a，擠壓速率為1.39 mm/a。這與張四昌［7］使用地質構造、地震活動圖像和震源機制解等資料綜合分析得到的結論一致。該部位相對較小的擠壓速率說明，豐臺-野雞坨斷裂作為唐山塊陷與薊縣塊隆的主要邊界，其擠壓效應不甚明顯，其變形情況與區域形變場具有相對較好的一致性，該斷裂帶活動性相對較弱。

2)薊運河斷裂:該斷裂主要以左旋擠壓運動為主，其中左旋走滑速率為6.31 mm/a，擠壓速率為5.27 mm/a。該位置擠壓速率較大，說明唐山塊體在受到北東方向上的山海關塊體的推擠下，在薊運河斷裂處受到阻礙，并進行快速的應變積累。該斷裂帶在現今構造環境下具有相對較大的活動性，地震危險性相對較大。

圖4 地殼運動模擬結果和GPS觀測結果的誤差比對(誤差棒代表GPS觀測結果的誤差范圍)Fig.4 Comparision of crustal movement by GPS observation with FEM simulation

圖5 唐山塊體內斷裂構造運動速度計算結果Fig.5 Movement rates of the faults in Tangshan block

3)唐山斷裂:主要呈現右旋走滑和拉張狀態，其右旋走滑速率為2.88 mm/a，拉張速率為0.21 mm/a。這與郭慧、江娃利等人［8］在唐山地區一系列鉆探和孫家樓槽探得到的結果具有一致性。其微弱的拉張速率說明唐山塊陷在鎖閉狀態下，塊體內部的唐山斷裂帶所處的構造環境相對穩定，暫時不具備發生大震的可能性。

4)寧河-昌黎斷裂:該斷裂主要以右旋拉張狀態為主，其走滑速率為1.25 mm/a，拉張速率為2.26 mm/a，說明樂亭塊陷對唐山塊陷具有明顯的拉張牽引作用。灤縣-樂亭斷裂主要呈現為擠壓逆沖狀態，其擠壓速率為1.4 mm/a，右旋運動速率為0.8 mm/a，說明山海關塊隆對唐山塊陷具有明顯的擠壓逆沖作用，這與虢順民［1］的認識一致。

表3 唐山地區主要斷裂帶運動速率模擬結果及與地質結論比較Tab.3 Comparison of movement rates from FEM simulation with geological analysis

4 結論

本文利用1999～2007年GPS數據和二維接觸有限元方法，對唐山地區斷裂帶運動進行數值模擬，計算了唐山塊陷內部及其周邊塊體內相關斷裂帶的運動狀態，得到以下主要結論:

1)總體而言，唐山地區地殼形變場具有較為明顯的不均勻性和分界性，總體呈現出向西運動的速度場狀態。研究區自北向南，除個別點以外，速度矢量方向由正西方向逐漸轉變為南西方向，不均勻性的同時，呈現較為規律的扇形分布。唐山塊陷內部GPS觀測點運動速率明顯小于樂亭塊陷、山海關塊隆和薊縣塊隆等外圍地塊，說明唐山塊陷在四周深大斷裂的圍陷下，在1999～2007年間呈現為相對鎖閉平衡狀態。

2)唐山斷裂帶在現今構造環境下的活動性相對較弱，不具備發生大震的危險性。但是山海關南西西向的推擠運動導致唐山塊體與冀中塊陷結合處的薊縣斷裂呈現出快速的應變積累狀態，其斷層活動速率較高，具有發生較大規模地震事件的可能性。

1 虢順民，李義志.唐山地震區域構造背景和發展模式的討論［J］.地質科學，1977(10):305 －321.(Guo Shunmin，Li Yizhi.Discussion on the regional structural background and the seismogenic model of the Tangshan earthquake［J］.Scientia Geologica sinica，1977(10):305 －321)

2 鄧起東，張培震，冉勇康，等.中國活動構造基本特征［J］.中國科學:D輯:地球科學，2002(12)1 020－1 030.(Deng Qidong，Zhang Peizhen，Ran Yongkang，et al.The basic characteristics of active tectonic structure［J］.Science in China:Series D，2002(12):1 020 －1 030)

3 雷顯權，陳運平，趙炯洋.天山現今地殼變形的非連續接觸模型模擬［J］.中南大學學報:自然科學版，2011(9):2 754 －2 762.(Lei Xianquan，Chen Yunping，Zhao Jiongyang.Discontinuous contact model to simulate current crustal deformation in Tianshan Mountain［J］.Journal of Central South University:Science and Technology，2011，42(9):2 754－2 762)

4 陳顒，黃庭芳.巖石物理學［M］.合肥:中國科學技術大學出版社，2009.(Chen Yong，Huang Tingfang.Rock physics［M］.Hefei:Press of University of Science and Technology of China，2009)

5 Vermilye J M，Scholz C H.The process zone:a microstructural view of faultgrowth［J］.Journal of Geophysical Research，1998，103(B6):12 223 –12 237.

6 安關峰，高大釗.接觸面彈粘塑性本構關系研究［J］.土木工程學報，2001，34(1):88 － 91.(An Guanfeng，Gao Dazhao.Research on elastic-visco-plastic constitution of interfaces［J］.Civil Engineering Journal，2001，34(1):88 －91)

7 張四昌，刁桂琴.唐山地震序列的構造過程［J］.中國地震，1992(2):75 － 82.(Zhang Sichang，Diao Guiqin.The tectonic process of the Tangshan earthquake sequence［J］.Earthquake Research in China，1992，8(2):73 －80)

8 郭慧，江娃利，謝新生.鉆孔與探槽揭示1976年河北唐山Ms7.8地震發震構造晚第四紀強震活動［J］.中國科學:地球科學，2011(7):1 009 －1 028.(Guo Hui，Jiang Wali，Xie Xinsheng.Late-Quaternary strong earthquakes on the seismogenic fault of the 1976 Ms7.8 Tangshan earthquake，Hebei，as revealed by drilling and trenching［J］.Sci China Earth Sci，2011，41(7):1 009 －1 028)

SIMULATING CURRENT CRUSTAL DEFORMATION IN TANGSHAN AREA WITH DISCONTINUOUS CONTACT MODEL CONSTRAINED BY GPS

Meng Qingxiao1，2，3)
1)First Crust Monitoring and Application Center，CEA，Tianjin 300180
2)Key Leboratory of Earthquake Geodesy，Institute of Seismology，CEA，Wuhan 430071
3)Wuhan Base of Institute of Crustal Dynamics，CEA，Wuhan430071

The current crustal deformation of the Tangshan area was analyzed with FEA numerical method under constraint of GPS observation，using physical contact model to obtain the discontinuous motion of faults.The results of simulation show that the movement of the faults in Tangshan area is different from each others.The faults in Tangshan area are right lateral slipping for Tangshan Block is pushed by Shanhaiguan plate.Jiyunhe fault zone has accumulated higher stain，therefore seismic risk there is higher.

Tangshan area;GPS;FEM;contact model;fault zone.

P227

A

1671-5942(2014)05-0040-05

2013-09-12

地震科技星火基金青年項目(XH14067Y);中國地震局地震研究所所長基金項目(IS201266113);中國地震局監測預報司震情跟蹤專項(2014020212)。

孟慶筱，男，1986年生，助理工程師，主要從事工程地震、水庫誘發地震、數值模擬和構造形變分析等研究。E－mail:mqxgucas@163.com。