李立功，段虎榮，牛麗娟

（1.陜西鐵路工程職業(yè)技術(shù)學(xué)院, 陜西 渭南 714000；2.西安科技大學(xué) 測繪科學(xué)與技術(shù)學(xué)院，陜西 西安 710054）

基于旋轉(zhuǎn)模型的青藏高原重力變化研究

李立功1，段虎榮2，牛麗娟1

（1.陜西鐵路工程職業(yè)技術(shù)學(xué)院, 陜西 渭南 714000；2.西安科技大學(xué) 測繪科學(xué)與技術(shù)學(xué)院，陜西 西安 710054）

用數(shù)值模擬方法研究了地殼旋轉(zhuǎn)運(yùn)動與空間重力變化的關(guān)系，得出地殼劃分的長方體的旋轉(zhuǎn)運(yùn)動引起的重力變化呈中心對稱分布。青藏高原地區(qū)的研究資料表明，該區(qū)域繞青藏高原東構(gòu)造存在順時針旋轉(zhuǎn)運(yùn)動。文中結(jié)合該地區(qū)的數(shù)字高程以及GPS觀測速率數(shù)據(jù)，計算了青藏高原塊體旋轉(zhuǎn)運(yùn)動引起的重力變化，并與GRACE衛(wèi)星重力觀測結(jié)果進(jìn)行了對比分析。

構(gòu)造；衛(wèi)星重力；旋轉(zhuǎn)運(yùn)動；重力變化

青藏高原是世界上海拔最高、最年輕、動力學(xué)環(huán)境最復(fù)雜的高原，是地學(xué)研究的熱點(diǎn)。衛(wèi)星重力探測技術(shù)以其快速、高效、幾乎覆蓋全球并能給出大尺度的重力變化等特點(diǎn)，大大改善了人們對地球重力場的認(rèn)識[1-4]。與傳統(tǒng)的重力觀測理論不同，現(xiàn)有的重力變化理論模型大都不太適合解釋衛(wèi)星重力觀測到的結(jié)果[5]。為此，本文建立了地殼旋轉(zhuǎn)模型，嘗試通過模型的旋轉(zhuǎn)運(yùn)動來模擬地殼運(yùn)動，進(jìn)而對引起的重力變化進(jìn)行分析解釋。

1 旋轉(zhuǎn)模型及計算方法

1.1 地殼旋轉(zhuǎn)產(chǎn)生的重力變化

由萬有引力定律可知，在重力場中，對于某個相對固定的觀測點(diǎn)（ε，η，ζ），點(diǎn)源物質(zhì)位移時產(chǎn)生的重力變化為[6]：

式中，G表示地球引力常數(shù)；m表示點(diǎn)源物質(zhì)的質(zhì)量；r1、r2分別表示位移前后觀測點(diǎn)到點(diǎn)源物質(zhì)的距離。

在進(jìn)行地質(zhì)體實際研究時,可將地質(zhì)體分成若干個規(guī)則的長方體[5]，每個長方體可以理解為由n個點(diǎn)源物質(zhì)組成。根據(jù)式（1），對第i個長方體運(yùn)動進(jìn)行積分，就可以得到該長方體發(fā)生位移時在空間某處產(chǎn)生的重力變化[6]。

式中，ρ表示介質(zhì)密度；V表示長方體的體積。

根據(jù)式（2），當(dāng)直立長方體做旋轉(zhuǎn)運(yùn)動時，為便于數(shù)據(jù)處理計算，不妨在建模時，將直立長方體的2個底面建立為正方形，設(shè)其底邊長為2l，棱長為2h。此時，在相應(yīng)空間點(diǎn)處的重力變化為：式中，θ為長方體旋轉(zhuǎn)運(yùn)動所轉(zhuǎn)過的角度，θ ∈[0,2π]。設(shè)（XC，YC，ZC）為長方體中心點(diǎn)坐標(biāo)，則X1，X2∈ [XC-l，XC+l]，Y1，Y2∈ [YC-l，YC+l]，Z1，Z2∈ [ZC-h，ZC+h]。

由重力的可加性知，所有長方體的旋轉(zhuǎn)運(yùn)動引起空間點(diǎn)（ε，η，ζ）處的重力變化為：

1.2 GRACE衛(wèi)星觀測數(shù)據(jù)計算重力變化[5,7]

GRACE衛(wèi)星采用高-低和低-低SST兩種模式進(jìn)行探測，可以精確測定地球重力場的中、長波部分的變化。在正常橢球情況下,對地面某點(diǎn)（r，θ，λ）處的地球重力位球函數(shù)r求導(dǎo)，可得：式中，Plm (c osθ)為完全規(guī)格化l階m級締合Legendre多項式；fM為地球引力常數(shù)為完全規(guī)格化的球函數(shù)系數(shù)；R為地球半徑；r為地面觀測點(diǎn)到地球質(zhì)心的距離；λ、θ為觀測點(diǎn)的經(jīng)度和緯度。

GARCE衛(wèi)星數(shù)據(jù)觀測得到的重力場球諧系數(shù)可以表示為重力場長期的靜態(tài)平均值、重力場的長期變化與固體潮、海潮、極潮的影響及大氣、非潮汐變化影響的組合。而對某個區(qū)域的某一時間跨度而言，重力場系數(shù)可以看成背景場與其改正之和[8]。

將GRACE衛(wèi)星每年獲得的重力場與背景重力場求差，就可得到年重力場的變化值：

即

其中，ΔC、ΔS 表示月重力場和背景重力場模型的球諧系數(shù)之差。為抑制GRACE衛(wèi)星重力模型球諧系數(shù)高階部分的噪聲及截斷誤差的影響，加入各項同性高斯平滑函數(shù)Wl即

2 數(shù)值模擬地殼旋轉(zhuǎn)位移引起的重力場變化

進(jìn)行數(shù)據(jù)模擬的單個長方體參數(shù)為：X軸方向（X1=0 km，X2=20 km），Y軸方向（Y1=0 km， Y2=20 km），Z軸方向（Z1=0 km, Z2=30 km）。計算點(diǎn)分布在ζ=0．0的平面上，共100個，點(diǎn)間隔10 km，運(yùn)動前后介質(zhì)密度ρ=2.67×103kg/m3。表1為長方體旋轉(zhuǎn)參數(shù)表。

表1 長方體旋轉(zhuǎn)參數(shù)表

圖1 旋轉(zhuǎn)運(yùn)動引起的重力異常等值線圖/10-8ms-2

3 青藏高原地區(qū)重力場變化

圖2中紅色線框范圍為本文研究區(qū)域。中國西部地區(qū)地殼運(yùn)動的位移場整體上具有由西南向北東方向運(yùn)動的特點(diǎn)，而且青藏高原地區(qū)繞喜馬拉雅東構(gòu)造具有明顯的順時針轉(zhuǎn)動模式。

圖2 青藏高原地區(qū)GPS運(yùn)動速度場

3.1 GRACE衛(wèi)星監(jiān)測青藏高原地區(qū)重力變化

受時間和空間分辨率的影響，GRACE衛(wèi)星在研究局部地區(qū)的重力變化及其原因上，效果不顯著。但是可以從大尺度上反映地球重力場的變化。又考慮到間隔時間長的旋轉(zhuǎn)運(yùn)動應(yīng)更明顯，根據(jù)式（6）和式（8），以青藏高原地區(qū)東經(jīng)80°～105°、北緯20°～35°地區(qū)GRACE GX-OG-_2-GSM 2005年12個月的重力場模型平均值作為背景重力場模型，與2009年12個月的重力場模型的平均值求差，計算2009年相對于2005年青藏高原地區(qū)的重力變化。模型最大階次為120階，高斯平滑半徑取390 km，分辨率為1°×1°[10]，如圖3所示。

圖3 2005～2009年青藏高原地區(qū)120階重力場變化/10-8ms-2

從圖中我們可以看到昆明以西地區(qū)是重力場變化值最大的地方，變化值為正值；而在喜馬拉雅山脈的南部以及青藏高原地區(qū)和川滇區(qū)域的交界部分，重力異常值為負(fù)值。

將上述結(jié)果與該地區(qū)的GPS觀測速度場（見圖2）進(jìn)行對比分析可以發(fā)現(xiàn)，因為青藏高原主體地殼向北運(yùn)動,而且呈現(xiàn)由南到北逐漸減小的趨勢，物質(zhì)的流出較為明顯, 所以該地區(qū)的重力變化值為負(fù)；而青藏高原地殼存在繞喜馬拉雅東構(gòu)造順時針旋轉(zhuǎn)運(yùn)動，受到四川地塊的阻擋作用，所以昆明以西地區(qū)的重力變化值為正[11]。

3.2 青藏高原地區(qū)數(shù)值模擬的重力變化

為方便計算，本文只考慮塊體作旋轉(zhuǎn)運(yùn)動的情況，建立如圖4所示坐標(biāo)系，把圖2所示研究區(qū)域分成15行25列，1°×1°大小的長方體。計算面為xoy平面（Z=0），即每個計算點(diǎn)都位于相應(yīng)長方體在xoy平面投影的中心位置。對于每個長方體塊體，地殼厚度統(tǒng)一取30 km，h表示塊體地形高度，hmax表示塊體地形的最大高度。

圖4 塊體地形高度

考慮到該地區(qū)的地形、氣候等條件的影響，從http://www．ngdc．noaa．gov/網(wǎng)站獲取青藏高原研究地區(qū)的地形數(shù)據(jù)，結(jié)合圖3中紅色線框所示范圍內(nèi)GPS觀測站的運(yùn)動速率[9]及《中國大陸現(xiàn)今地殼運(yùn)動》附錄二中記錄的中國大陸GPS測站的速度場來確定塊體運(yùn)動[12]。對于GPS觀測點(diǎn)分布不均勻的塊體，采用以計算點(diǎn)為中心，以1°為半徑搜索，按距離平方的倒數(shù)定權(quán)的方法進(jìn)行加權(quán)平均內(nèi)插得到。塊體的旋轉(zhuǎn)運(yùn)動量為塊體內(nèi)所有旋轉(zhuǎn)運(yùn)動量的平均值，由公式（3）和（4）可得青藏高原地區(qū)重力變化如圖5。

將旋轉(zhuǎn)數(shù)據(jù)計算的重力變化與應(yīng)用GRACE衛(wèi)星數(shù)據(jù)計算的重力變化進(jìn)行綜合對比分析可知，兩者的圖形形態(tài)和數(shù)值都有較好的一致性。青藏高原地區(qū)的重力變化為負(fù)，昆明以西地區(qū)的重力變化值為正。但在部分地區(qū)也存在差異，這些差異需進(jìn)一步研究。

圖5 由旋轉(zhuǎn)運(yùn)動引起的青藏高原地區(qū)重力變化/10-8ms-2

[1] Trung N N,Lee S M,Que B C．Satellite gravity Anomalies and Their Correlation with the Major Tectonic Features in the South China Sea [J]． Gortdwand Research,2004,7(2):407-424

[2] Han S,Shum C,Bevis M,etal．Crustal Dilatation Observed by GRACE after the 2004 Sumatra-Andam an Earthquake[J]．Science,2006,313(5 787):658-662

[3] 周旭華,許厚澤,吳斌,等．用GRACE衛(wèi)星跟蹤數(shù)據(jù)反演地球重力場[J]．地球物理學(xué)報,2006,49(3):718-723

[4] 鄭偉,許厚澤,鐘敏,等．國際重力衛(wèi)星研究進(jìn)展和我國將來衛(wèi)星重力測量計劃[J]．測繪科學(xué),2010,35(1):5-9

[5] 段虎榮,張永志,徐海軍,等．中國西部地殼垂直運(yùn)動引起的重力場空間變化特征[J]．大地測量與地球動力學(xué), 2011,31(3):25-28

[6] 段虎榮,張永志,徐海軍,等．中國西部地殼水平運(yùn)動引起的重力場空間變化特征[J]．大地測量與地球動力學(xué),2011,31(1):24-28

[7] 寧津生．衛(wèi)星重力探測技術(shù)與地球重力場研究[J]．大地測量與地球動力學(xué),2002,22(1):1-5

[8] 李立功．基于旋轉(zhuǎn)模型的青藏高原重力變化研究[D]．西安:長安大學(xué),2011

[9] 王琪,張培震,馬宗晉．中國大陸現(xiàn)今構(gòu)造變形GPS觀測數(shù)據(jù)與速度場[J]．地學(xué)前沿,2002,9(2):415-429

[10] 段虎榮,張永志,劉鋒,等．利用衛(wèi)星重力數(shù)據(jù)研究中國及鄰域地殼厚度[J]．地球物理學(xué)進(jìn)展，2010,25(2):494-499

[11] 曹建玲,石耀霖,張懷,等．青藏高原GPS位移繞喜馬拉雅東構(gòu)造結(jié)順時針旋轉(zhuǎn)成因的數(shù)值模擬[J]．科學(xué)通報,2009,54(2):224-234

[12] 賴錫安,黃立人,徐菊生．中國大陸現(xiàn)今地殼運(yùn)動[M]．北京:地震出版社,2004

P223.0

B

1672-4623（2014）04-0037-03

10.11709/j.issn.1672-4623.2014.04.012

李立功，碩士，研究方向為衛(wèi)星重力。

2013-10-09。

項目來源：陜西省教育廳自然科學(xué)基金資助項目（12JK0798）；陜西省科學(xué)技術(shù)廳自然科學(xué)基礎(chǔ)研究計劃資助項目（2011JM5005）；陜西鐵路工程職業(yè)技術(shù)學(xué)院科研基金資助項目（2011-40）。