數(shù)值模擬華北克拉通巖石圈減薄的一種可能機(jī)制——下地殼榴輝巖重力失穩(wěn)引起的拆沉

喬彥超，郭子祺，石耀霖

1 中國(guó)科學(xué)院遙感應(yīng)用研究所，北京 100101

2 中國(guó)科學(xué)院研究生院，北京 100039

1 引 言

大陸巖石圈一直被認(rèn)為是相對(duì)穩(wěn)定的構(gòu)造單元，然而并不是所有的前寒武紀(jì)克拉通從它們形成以后都是處于穩(wěn)定狀態(tài).近期的研究發(fā)現(xiàn)，如北美的Wyoming克拉通自中元古代以來(lái)發(fā)生了減薄和破壞［1－2］，南非的 Kapvaal克拉通在3.1Ga 時(shí)由于巖石圈地幔的拆沉作用而被改造［3］，Seber等研究發(fā)現(xiàn)地中海Alboran海盆及其周圍摩洛哥北部Rif造山帶和西班牙南部Betic造山帶（兩者均屬阿爾卑斯造山期產(chǎn)物）巖石圈地?，F(xiàn)今正在發(fā)生拆沉作用［4］，很多結(jié)果也表明［5－18］華北克拉通在古生代再活化而發(fā)生減薄，克拉通巖石圈的改造和破壞已經(jīng)成為大陸動(dòng)力學(xué)研究的一個(gè)熱點(diǎn).也正是如此，華北克拉通巖石圈破壞已成為繼青藏高原和大別—蘇魯超高壓變質(zhì)帶之后，我國(guó)又一在國(guó)際地球動(dòng)力學(xué)界引起廣泛關(guān)注的重大科學(xué)問(wèn)題，已成為國(guó)際地球科學(xué)前沿研究領(lǐng)域.

關(guān)于華北克拉通巖石圈減薄機(jī)制還是有比較多的爭(zhēng)議［6－8，12，17－20］，高山等［17－18］通過(guò)研究在華北克拉通發(fā)現(xiàn)高鎂安山巖、英安巖、和埃達(dá)克巖，并認(rèn)為這些巖石是部分熔融的榴輝巖和地幔巖石相互反應(yīng)生成，是華北克拉通拆沉的主要證據(jù).拆沉作用泛指由于重力的不穩(wěn)定性導(dǎo)致巖石圈地幔、大陸下地殼或大洋地殼沉入下伏軟流圈或地幔的過(guò)程.俯沖大洋巖石圈中榴輝巖的拆沉再循環(huán)是板塊構(gòu)造的直接產(chǎn)物［21－23］.

拆沉理論認(rèn)為巖石圈加厚是大陸巖石圈和下地殼拆沉的前奏.巖石圈加厚促使基性下地殼轉(zhuǎn)變?yōu)榱褫x巖［17－18］；榴輝巖的密度比巖石圈地幔的密度高200～400kg／m3［24－26］.這種差異將造成重力上的不穩(wěn)定性，使榴輝巖拆沉再循環(huán)進(jìn)入地幔中［27－32］.

本文研究通過(guò)數(shù)值計(jì)算的方法，結(jié)合華北克拉通已有的地質(zhì)成果，模擬華北克拉通巖石圈拆沉機(jī)制，通過(guò)研究討論能夠發(fā)生拆沉需要榴輝巖的厚度，不同規(guī)模榴輝巖造成拆沉的不同模式，以及對(duì)比熱侵蝕和拆沉所造成的地表不同反應(yīng)，來(lái)揭示華北克拉通巖石圈不同減薄機(jī)制的動(dòng)力所對(duì)應(yīng)的不同地質(zhì)構(gòu)造特征.

2 模 型

2.1 數(shù)值計(jì)算的基本原理

我們?cè)鴮?duì)計(jì)算原理、使用公式以及所使用程序的驗(yàn)證進(jìn)行過(guò)論述［23］，本文不再詳述.研究基于粘滯系數(shù)、密度與溫度相關(guān)的二維上地幔小尺度地幔對(duì)流［33－40］，其關(guān)系式為

其中η0是等效粘滯系數(shù)，b是經(jīng)驗(yàn)系數(shù)，Ttop和Tbottom分別是計(jì)算區(qū)域上下邊界的溫度，ρ0是等效密度，α是膨脹系數(shù).

我們使用了有限差分方法，并結(jié)合MIC技術(shù)［34－36］求解連續(xù)性方程，動(dòng)量守恒方程和能量守恒方程.計(jì)算中使用無(wú)限大Prandtl數(shù)假定，這樣忽略動(dòng)量方程中的慣性項(xiàng)，假定流體是不可壓縮的，并且引入 Boussinesq假設(shè)［41，42］.

計(jì)算的二維模型是長(zhǎng)700km，寬700km的方盒子，網(wǎng)格為71×71，共有400×400個(gè)示蹤點(diǎn).邊界條件：所有的邊都為自由邊界條件；左右兩個(gè)邊界是絕熱邊界條件 （?T／?x ＝0）；頂邊界溫度為Ttop＝273K；底邊界溫度為Tbottom＝1773K.

克拉通拆沉過(guò)程中的第一步是巖石圈加厚玄武巖相變?yōu)槊芏缺容^大的榴輝巖.計(jì)算中我們使用由Wylllie［43］提出，Anderson［44］修改的玄武巖－榴輝巖相變圖，如圖1.因此大于50km的地殼被認(rèn)為是榴輝巖質(zhì)的，小于25km是玄武質(zhì)，深度在25～50km之間且溫度達(dá)到Teclo（z）（公式3，圖1）的玄武巖轉(zhuǎn)變?yōu)榱褫x巖.

圖1 巖石圈加厚玄武巖相變?yōu)榱褫x巖相變圖Fig.1 Phase transition of basalt to eclogite when lithospheres thicken

當(dāng)溫度達(dá)到固相線時(shí)巖石開始融化，計(jì)算溫度，選擇 McKenzie［45］公式（圖1）：

本文使用的是溫度場(chǎng)定義的巖石圈厚度.巖石圈加厚前的地溫曲線為圖1中黑線，上地殼地溫梯度為25℃／km，則20km處為500℃；巖石圈100km厚，溫度為1350℃［46］；最底邊200km 處溫度為1400℃，如圖1中Told.地殼厚度為35km，即圖中原來(lái)Moho的深度.假定巖石圈快速擠壓增厚一倍到200km深，那么新的地溫曲線為圖1中綠色曲線，上地殼變?yōu)?0km，溫度為500℃；巖石圈變?yōu)?00km，溫度為1350℃，如圖1中Tnew.Tnew與Telo的交點(diǎn)深度約為34km，加厚的Moho深度為70km，那么從34km到70km都為相變區(qū)域，我們首先考慮極端情形，假定這一區(qū)域都變?yōu)榱褫x巖.

2.2 模型初始溫度場(chǎng)和物質(zhì)組成

我們計(jì)算區(qū)域?yàn)?00km×700km，溫度場(chǎng)主要參考Tnew的值，在底邊界700km處溫度為1500℃，如圖2.關(guān)于物質(zhì)組成，我們考慮極端的情形，數(shù)值模擬中我們往往對(duì)現(xiàn)實(shí)的模型或者資料進(jìn)行很多簡(jiǎn)化或者合理的假設(shè)，從而專注于我們研究的問(wèn)題，并且數(shù)值模擬的一大優(yōu)勢(shì)是，在保證正確的前提下，我們可以進(jìn)行，不同參數(shù)，不同條件，不同模型的討論.使我們得到研究問(wèn)題的全面資料.本研究中我們假定初始地殼厚度為70km，榴輝巖厚度為34～70km，首先計(jì)算橫向榴輝巖規(guī)模200～500km，如圖3.我們將討論橫向榴輝巖不同規(guī)模，如300～400km，200～500km，100～600km，即水平尺度分別為100km，200km，300km時(shí)所造成不同的拆沉模式.

圖2 計(jì)算的初始溫度場(chǎng)Fig.2 Initial thermal state

圖3 初始物質(zhì)組成Fig.3 Initial density

2.3 參數(shù)

計(jì)算中使用的參數(shù)如表1.計(jì)算中關(guān)鍵的參數(shù)是不同組分密度的選取，榴輝巖的密度取3500kg／m3［47］，地殼的密度取2800kg／m3，地幔的密度取3340kg／m3［48］.

表1 計(jì)算中用到的參數(shù)Table 1 Parameters used in computation

3 計(jì)算結(jié)果

榴輝巖的厚度為34～70km，我們橫向上采用三種不同規(guī)模的榴輝巖分別為200～500km，300～400km，100～600km，討論了不同規(guī)模的榴輝巖造成的不同拆沉結(jié)果.其溫度場(chǎng)等值線不同時(shí)刻結(jié)果如圖4.

圖4中第一橫排為中等規(guī)模榴輝巖（200～500km）情況下，不同時(shí)刻溫度場(chǎng)拆沉的等溫線結(jié)果，依次第二橫排為較小規(guī)模榴輝巖（300～400km），第三橫排為較大規(guī)模榴輝巖（100～600km）.通過(guò)對(duì)比這三種情況，我們可以發(fā)現(xiàn)，榴輝巖的規(guī)模越大拆沉的越劇烈，最終能發(fā)生拆沉的區(qū)域也越大.為了更好的觀察拆沉的過(guò)程，同時(shí)我們給出了不同模型密度的隨時(shí)間的變化過(guò)程如圖5.

圖5的橫排次序如圖4.從圖5中可見，當(dāng)榴輝巖規(guī)模較大時(shí)，拆沉開始階段的過(guò)程更傾向于雙管道向下拆沉，模型較短時(shí)為單管道拆沉的模式.為了得到更精確的定量數(shù)據(jù)，我們給出了三種模型的榴輝巖最終時(shí)刻的溫度彩圖疊加速度矢量圖.其結(jié)果如圖6—8.

圖6中三條黑色線從上到下分別為500，1000，1350℃時(shí)的溫度等值線圖，綠色線為拆沉前的巖石圈底邊界，即200km處.因此從綠色線到最下面的黑線之間的區(qū)域?yàn)槿A北克拉通巖石圈能夠發(fā)生拆沉的區(qū)域.通過(guò)計(jì)算我們能夠知道減薄后的巖石圈最大厚度為64km，發(fā)生拆沉的最大厚度為136km.從流體場(chǎng)可以看出最終時(shí)刻發(fā)育兩個(gè)比較完整的對(duì)流環(huán).整個(gè)流體場(chǎng)比較穩(wěn)定.

圖7中的線代表的意義如圖6.通過(guò)計(jì)算我們能夠知道減薄后的巖石圈最大厚度為95km，發(fā)生拆沉的最大厚度為105km.可以看到最終時(shí)刻1350℃等溫線以上兩端局部有小的對(duì)流環(huán)，中心區(qū)域則是比較穩(wěn)定，這樣兩端可能會(huì)熱流更高，對(duì)應(yīng)著比較強(qiáng)烈，時(shí)間比較長(zhǎng)的巖漿活動(dòng).同華北克拉通巖漿分布地質(zhì)資料［51］相比也比較符合.

圖8中的線代表的意義如圖6.通過(guò)計(jì)算我們能夠知道減薄后的巖石圈最大厚度為108km，發(fā)生拆沉的最大厚度為92km.可以看到最終時(shí)刻流體場(chǎng)發(fā)育多管道向下拆沉.

圖4 三種模型不同時(shí)刻溫度場(chǎng)等值線圖Fig.4 Temperature contours of three models at different times

圖5 三種模型不同時(shí)刻密度場(chǎng)變化圖Fig.5 Density variations of three models at different times

圖6 較長(zhǎng)模型時(shí)最終時(shí)刻拆沉的結(jié)果Fig.6 Delamination result of the big model at last

圖7 模型為中等時(shí)最終時(shí)刻拆沉的結(jié)果Fig.7 Delamination result of the middle model at last

圖8 模型為最小時(shí)最終時(shí)刻拆沉的結(jié)果Fig.8 Delamination result of the small model at last

4 拆沉減薄與熱侵蝕減薄的對(duì)比討論

通過(guò)對(duì)華北克拉通巖石圈熱對(duì)流減薄機(jī)制的計(jì)算［32］及本文對(duì)華北克拉通巖石圈拆沉減薄機(jī)制的研究對(duì)比，我們可以知道在力學(xué)機(jī)制上這兩種機(jī)制都是可能發(fā)生的，即華北克拉通巖石圈由于熱擾動(dòng)引起對(duì)流減薄和由于榴輝巖密度不穩(wěn)定引起對(duì)流減薄都是可能發(fā)生的.但是通過(guò)細(xì)致的對(duì)比如圖9，我們可以發(fā)現(xiàn)不同的減薄機(jī)制的結(jié)果是不一致的，對(duì)應(yīng)的地表地質(zhì)特征也完全不同.

圖9 熱對(duì)流減薄和拆沉減薄結(jié)果對(duì)比圖Fig.9 Comparison of heat convection thinning and delamination

從圖9中左圖為熱對(duì)流減薄的結(jié)果，右圖是拆沉減薄的結(jié)果，從圖中藍(lán)色箭頭標(biāo)志我們可以發(fā)現(xiàn)，熱對(duì)流減薄，中心區(qū)域變熱變薄，中心區(qū)域是一個(gè)伸展區(qū)域，發(fā)展一個(gè)大的變質(zhì)核雜巖區(qū)，巖漿的發(fā)展時(shí)序是從中心向兩側(cè)對(duì)稱分布，并且向兩側(cè)減弱；熱拆沉減薄結(jié)果中，中心區(qū)域受到擠壓加厚，反而變冷，但是我們知道，對(duì)于地質(zhì)歷史上發(fā)生過(guò)的加熱事件通過(guò)地質(zhì)手段，如探針等方法是可以得到信息的，對(duì)于變冷的過(guò)程我們是沒(méi)有辦法得到信息的.在中心區(qū)域的兩側(cè)發(fā)生兩個(gè)巖漿活動(dòng)中心，巖漿的發(fā)展時(shí)序是從兩側(cè)向中心減弱發(fā)展.

通過(guò)上面討論，何種減薄機(jī)制在這一過(guò)程中起主導(dǎo)作用，要更多的分析地表的地質(zhì)資料，來(lái)找到足夠的證據(jù).

5 結(jié) 論

本文研究表明，中生代華北拉通巖石圈受擠壓，玄武巖相變?yōu)槊芏容^大的榴輝巖，榴輝巖重力失穩(wěn)能夠引起下地殼的拆沉，造成巖石圈的大規(guī)模減薄.對(duì)比不同規(guī)模的榴輝巖減薄方式發(fā)現(xiàn)，當(dāng)榴輝巖較大規(guī)模時(shí)更可能發(fā)生雙管道拆沉，而規(guī)模小時(shí)發(fā)生單管道拆沉.榴輝巖模型最大時(shí)候（100～600km）減薄后的巖石圈最大厚度為64km，發(fā)生拆沉的最大厚度為136km；規(guī)模中等時(shí)候（200～500km）巖石圈最大厚度為95km，發(fā)生拆沉的最大厚度為105km；規(guī)模比較小時(shí)候（300～400km）巖石圈最大厚度為108km，發(fā)生拆沉的最大厚度為92km.

通過(guò)對(duì)比中生代華北克拉通巖石圈熱減薄和拆沉減薄機(jī)制可知，從力學(xué)機(jī)制上都能造成巖石圈減薄，但他們的結(jié)果對(duì)應(yīng)著不同的地質(zhì)特征.對(duì)比我們的計(jì)算，熱侵蝕最多能減薄100km［32］，拆沉最多能夠減薄136km，但是很多資料都表明華北東部最多減薄了120km以上，所以單從減薄量上來(lái)看拆沉減薄機(jī)制更符合實(shí)際地質(zhì)資料.根據(jù)華北克拉通巖漿分布規(guī)律［7，49－50］我們可以發(fā)現(xiàn)有東北部和南部?jī)蓚€(gè)比較大的中心，這與我們文中討論的中等規(guī)模的榴輝巖拆沉減薄機(jī)制能夠在中心兩側(cè)發(fā)展兩個(gè)巖漿活動(dòng)中心比較一致.同時(shí)如果我們單純考慮熱侵蝕減薄機(jī)制，它所對(duì)應(yīng)的底部熱擾動(dòng)太過(guò)于偶然性，為什么不在其他地方，而只在華北克拉通，相比較于熱侵蝕，由于拆沉理論得到了很多地質(zhì)學(xué)的證據(jù)，如高山等的地球化學(xué)證據(jù)證明華北發(fā)生拆沉進(jìn)入地幔的榴輝巖存在，所以我們更傾向于榴輝巖重力失穩(wěn)引發(fā)拆沉是主導(dǎo)華北克拉通巖石圈減薄的機(jī)制.

（References）

［1］Sleep N H.Evolution of the continental lithosphere.Annual Review of Earth and Planetary Sciences，2005，33：369－393.

［2］Carlson R W，Pearson D G，Physical，chemical and chronological characteristics of continental mantle.Reviews of Geophysics，2005，43：RG1001，doi：10.1029／2004RG000156.

［3］Moser D E，F(xiàn)lowers R M，Hart R J.2001.Birth of the Kaapvaal tectosphere 3.08billion years ago.Science，291：465－468.

［4］Seber D，Barazangi M，lbebrahim A.Geophysical evidence for lithospheric delamination beneath the Alboran Sea and Rif－Betic mountains.Nature，1996，379：785－790.

［5］Chen L，Zheng T，Xu W.A thinned lithospheric image of the Tanlu Fault Zone，eastern China：Constructed from wave equation based receiver function migration.Journal of Geophysical Research，2006，111：B09312.

［6］高山，章軍鋒，許文良，等.拆沉作用與華北克拉通破壞.科學(xué)通報(bào)，2009，54：1962－1973.Gao S，Zhang J F，Xu W L，et al.Delamination and destruction of the North China Craton.Chinese Sci Bull，2009，54，doi：10.1007／s11434－009－0395－9.

［7］吳福元，徐義剛，高山，等.華北巖石圈減薄與克拉通破壞研究的主要學(xué)術(shù)爭(zhēng)論.巖石學(xué)報(bào)，2008，24：1145－1174.Wu F Y，Xu，Y G，Gao S，et al.Controversy over studies of the lithospheric thinning and craton destruction of North China.Acta Petrol Sin （in Chinese），2008，24：1145－1174.

［8］郝曉光，方劍，劉根友等.巖石圈拆沉還是軟流圈上涌？——正常密度理論與華北克拉通破壞重力學(xué)機(jī)制研究.地球物理學(xué)進(jìn)展，2011，26（4）：1137－1140，doi：10.3969／j.issn.1004－2903.2011.04.001.Hao X G ，F(xiàn)ang J，Liu G Y ，et al.North China craton destruction，delamination or thermal erosion？Chinese J.Geophys.（in Chinese），2011，26（4）：1137－1140，doi：10.3969／j.issn.1004－2903.2011.04.001 .

［9］許衛(wèi)衛(wèi)，鄭天愉，趙亮.華北地區(qū)410km間斷面和660km間斷面結(jié)構(gòu)—克拉通活化的地幔動(dòng)力學(xué)狀態(tài)探測(cè).2011，中國(guó)科學(xué)：地球科學(xué)，41：678–685.Xu W W，Zheng T Y，Zhao L.Mantle dynamics of the reactivating North China Craton：Constraints from the topographies of the 410km and 660km discontinuities.Sci China Earth Sci，2011，54：881－887，doi：10.1007／s11430－010－4163－0.

［10］張旗.從"巖石圈減薄"到"克拉通破壞".地球物理學(xué)進(jìn)展，2011，26（6）：2262－2269，doi：10.3969／j.issn.1004－2903.2011.06.047.Zhang Q.The evolution from"lithospheric thinning"to"cratonic destruction" .Chinese J.Geophys.（in Chinese），2011，26（6）：2262－2269，doi：10.3969／j.issn.1004－2903.2011.06.047.

［11］Fan W M，Zhang H F，Baker J.On and off the north China craton：Where is the Archaean keel？Journal of Petrology，2000，41：933－950.

［12］Xu Y G. Thermo－tectonic destruction of the Archean lithospheric keel beneath the Sino－Korean Craton in China：Evidence，timing and mechanism.Physics and Chemisty of the Earth，2001，26：747－757.

［13］Liu Y S，Gao S.Melt－peridotite interactions：Links between garnet pyroxenite and high－Mg signature of continental crust.Earth and Planetary Science Letters，2005，234：39－57.

［14］Zheng J P，Griffin W L，O′Reilly S Y.Late Mesozoic－Eocene mantle replacement beneath the eastern North China Craton：Evidence from the Paleozoic and Cenozoic peridotite xenoliths.International Geology Review，2005，47：457－472.

［15］Zheng J P，Sun M，Zhou M F.Trace elemental and PGE geochemical constraints of Mesozoic and Cenozoic peridotitic xenoliths on lithospheric evolution of the North China Craton.Geochimica et Cosmochimica Acta，2005，69：3401－3418.

［16］Zheng J P，Griffin W L，O′Reilly S Y.Mineral chemistry of garnet peridotites from Paleozoic，Mesozoic and Cenozoic lithosphere：Constraints on mantle evolution beneath eastern China.Journal of Petrology，2006，47：2233－2256.

［17］Gao S，Rudnick R L，Carlson R W.Re－Os evidence for replacement of ancient mantle lithosphere beneath the North China craton.Earth and Planetary Science Letters，2002，198：307－322.

［18］Gao S，Rudnick R L，Yuan H L.Recycling lower continental crust in the North China craton.Nature，2004，432：892－897.

［19］Zhang H F，Sun M，Zhou X H.Mesozoic lithosphere destruction beneath the North China Craton：Evidence from major，trace element，and Sr－Nd－Pb isotope studies of Fangcheng basalts.Contributions to Mineralogy and Petrology，2002，144：241－253.

［20］Zhang H F，Sun M，Zhou X H.Secular evolution of the lithosphere beneath the eastern North China Craton：Evidence from Mesozoic basalts and high－Mg andesites.Geochimica et Cosmochimica Acta，2003，67：4373－4387.

［21］Sobolev A V，Hofmann A W，Sobolev S V，et al.An olivine－free mantle source of Hawaiian shield basalts.Nature，2005，434：590－597.

［22］Sobolev A V，Hofmann A W，Kuzmin D V，et al.The amount of recycled crust in sources of mantle－derived melts.Science，2007，316：412－417.

［23］郭慧麗，徐佩芬.地震層析成像在華北克拉通地區(qū)的研究進(jìn)展.地球物理學(xué)進(jìn)展，2011，26（5）：1557－1565，doi：10.3969／j.issn.1004－2903.2011.05.007.Guo H L，Xu P F.Progress of Seismic tomography applied in the North China Craton.Chinese J.Geophys.（in Chinese），2011，26（5）：1557－1565，doi：10.3969／j.issn.1004－2903.2011.05.007.

［24］Jull M，Kelemen.On the conditions for lower crustal convective instability.Journal of Geophysical Research，2001，106，6423－6446.

［25］Anderson D A.Large igneous provinces，delamination，and fertile mantle.Elements，2006，1：271－275.

［26］Anderson D A.Speculations on the nature and cause of mantle heterogeneity.Tectonophysics，2006，146：7－22.

［27］Lustrino M.How the delamination and detachment of lower crust can influence basaltic magmatism. Earth－Science Review，2005，72：21－38.

［28］Bédard J H. A catalytic delamination－driven model for coupled genesis of Archaean crust and sub－continental lithospheric mantle.Geochimica et Cosmochimica Acta，2006，70：1188－1214.

［29］Arndt N T，Goldstein S L.An open boundary between lower continental crust and mantle：Its role in crust formation and crustal recycling.Tectonophysics，1989，161：201－212.

［30］Kay R W，Kay S M.Creation and destruction of lower continental crust.Geologische Rundschau，1991，80：259－278.

［31］Rudnick R L.Making continental crust.Nature，1995.378：571－577.

［32］喬彥超，郭子祺，石耀霖.數(shù)值模擬華北克拉通巖石圈熱對(duì)流侵蝕減薄機(jī)制.中國(guó)科學(xué)：地球科學(xué)，2012，（待發(fā)表）.Qiao Y C，Guo Z Q，Shi Y L.Thermal Convection Thinning of the North China Craton：Numerical Simulation.Sci China Earth Sci.（in press）.

［33］Huang J S，Zhong S，van Hunen J.Controls on sublithospheric small－scale convection.Journal of Geophysical Research，2003，108（B8），2405.Doi：10.1029／2003JB002456.

［34］Gerya T V，Yuen D A.Characteristics－based marker－in－cell method with conservative finite－differences schemes for modeling geological flows with strongly variable transport properties.Physics of the Earth and Planetary Interiors，2003，140：295－320.

［35］Gerya T V， Maresch W V， Willner A P.Inherent gravitational instability of thickened continental crust with regionally developed low－timedium－pressure granulite facies metamorphism.Earth and Planetary Science Letters，2001，190：221－235.

［36］Gerya T V，Perchuk L L，Maresch W V.Thermal regime and gravitational instability of multi－layered continental crust：implications for the buoyant exhumation of high－grade metamorphic rocks.European Journal of Mineralogy，2002，14：687－699.

［37］葉正仁，王建.上地幔變黏度小尺度對(duì)流的數(shù)值研究.地球物理學(xué)報(bào)，2003，46（3）：335－339.Ye Z R，Wang J.A numerical research on the small scale convection with variable viscosity in the upper mantle.Chinese J.Geophys.（in Chinese），2003，46（3）：335－339.

［38］傅容珊，常筱華，黃建華，劉文忠.區(qū)域重力均衡異常和上地幔小尺度對(duì)流模型.地球物理學(xué)報(bào)，1994，37（S2）：249－258.Fu R S，Chang X H，Huang J H，et al.Regional gravity isostatic anomaly and small scale convection model in upper mantle.Chinese J.Geophys.（in Chinese），1994，37（S2）：249－258.

［39］熊熊，傅容珊，許厚澤等.增厚大陸巖石層熱邊界層對(duì)流剝離的數(shù)值模擬.地球物理學(xué)報(bào)，1998，41（S1）：33－40.Xiong X，F(xiàn)u R S，Xu H Z，et al.Numerical Simulation of the convective removal of the thickened thermal boundary layer of continental lothosphere.Chinese J.Geophys.（in Chinese），1998，41（S1）：33－40.

［40］傅容珊，黃建華，董樹謙.利用地震層析成象數(shù)據(jù)計(jì)算地幔對(duì)流新模型的探討.地球物理學(xué)報(bào)，2003，46（6）：772－778.Fu R S.A new mantle convection model constrained by seismic tomography data.Chinese J.Geophys.（in Chinese），2003，46（6）：772－778.

［41］Boussinesq J.Theorie analytique de la chaleur mise en harmonie avec la thermodynamique et avec la theorie mecanique de la lumiere.Gauthier－Villars Paris，1903，2：157－176.

［42］Tackley P. Effect of strongly temperature－dependant viscosity on time－dependent 3－dimensional model of mantle convection.Journal of Geophysical Research，1993，20（20）：2187－2190.

［43］Wyllie P.The Dynamic Earth.New York：Wiley，1971：416.

［44］Anderson D L. Tectonics and composition of Venus.Geophys.Res.Lett.，1980，1：101－102.

［45］McKenzie D.The generation and compaction of partially molten rock.Journal of Petrology，1984，25（3）：713－765.

［46］McKenzie D，Roberts J M，Weiss N O.Convection in the earth＇s mantle：Towards a numerical simulation.Journal of Fluid Mechanics，1974，6：465－538.

［47］Anderson D L.The upper mantle transition region：eclogite？Geophys.Res.Lett.，1979，6：433－435.

［48］Dupeyrat L.Thermal and chemical convection in planetary mantles.Journal of Geophysical Research，1995，100：497－520.

［49］吳福元，李獻(xiàn)華，楊進(jìn)輝等.花崗巖成因研究的若干問(wèn)題.巖石學(xué)報(bào)，2007，23（4）：1217－1238.Wu F Y，Li X H，Yang J H，et al.Discussions on the petrogenesis ofgranites.Acta Petrol.Sin.（in Chinese），2007，23（4）：1217－1238.

［50］嵇少丞，王茜，許志琴.華北克拉通破壞與巖石圈減薄.地質(zhì)學(xué)報(bào)，2008，82（2）：174－193.Ji S C，Wang Q，Xu Z Q.Break－up of the North China Craton through lithospheric thinning.Acta Geologica Sinica（in Chinese），2008，82（2）：174－193.致 謝 感謝鄭亮先生提供并維護(hù)并行機(jī)器.