造山带下岩石圈地幔小尺度减薄的动力学成因

地质与地球物理观测数据表明青藏高原、安第斯山、以及帕米尔等典型造山高原之下均有明显的岩石圈地幔小尺度/分段式减薄现象.这些小尺度岩石圈减薄难以用经典的拆沉或对流减薄理论来解释,一方面,拆沉预示大尺度岩石圈地幔的剥离过程,而对流减薄则在黏度相对低的地幔岩石圈中发生,其主要以小尺度的局部增厚触发并仅减薄地幔岩石圈的底部区域.另一方面,拆沉或对流减薄模型都预测造山带尺度的地幔岩石圈拆离,都假设造山带岩石圈横向均一,然而实际的造山带岩石圈往往由多个不同的地块构成,块体之间岩性、物性、流变结构可能大有差别,即横向不均一性.这些造山带岩石圈地幔的横向不均一性,能否有效解释观测到的局部小尺度减薄现象？为此,我们构建了一系列高精度动力学数值模型,系统模拟了碰撞造山过程中岩石圈地幔的形变和不稳定性.结果表明,在塑性屈服强度很低的情况下,横向不均一的造山带岩石圈有发生分段式/小尺度减薄的可能性;其主要机理是由位错蠕变与强塑性作用所导致的应变集中使得地块间及壳幔间耦合弱化,从而使得较弱地块的岩石圈地幔在增厚时由于重力不稳定性而产生局部剥离,进而诱发小尺度软流圈上涌.模拟结果可以良好地解释发生在青藏高原东北缘、安第斯中部高原、以及帕米尔高原之下岩石圈的局部小尺度/分段式减薄现象.     

榴辉岩下地壳拆沉的数值模拟

镁铁质的下地壳在一定的条件下可相变为榴辉岩并拆沉进入地幔,榴辉岩下地壳的拆沉可造成大陆地壳向长英质演化.地球化学研究显示,秦岭大别山地区在碰撞后期大约145-130Ma时还存在着加厚的榴辉岩下地壳;而到了130Ma时,加厚下地壳消失并产生了铁镁质和花岗质岩浆侵入.碰撞造山阶段,加厚的冷的地幔岩石圈的重力不稳定可造成岩石圈的拆沉,并携带下地壳一起进入地幔.但对碰撞造山构造运动后期的造山带或克拉通区域,地幔岩石圈基本处于稳定状态,下地壳榴辉岩化将对下地壳的拆沉起重要作用.本文采用二维有限元数值模拟的方法研究了岩石圈地幔重力稳定状态下,下地壳榴辉岩的拆沉过程以及拆沉时间.模型的主要参数包括下地壳榴辉岩的黏性(η_2)、密度(B_2)、总体规模(s)和岩石圈地幔的黏性(η_1).数值计算结果显示,由于下地壳榴辉岩的重力不稳定性,它能与其下的稳定岩石圈一起发生拆沉,拆沉过程的持续时间(t)与模型参数密切相关,大量计算分析显示t=7.3005×10~(-7)exp(0.6593η_2/η_1)η_1~(1.066)B_2~(-0.688)s~(0.1451η_2/η_1-0.9831).结果表明,岩石圈地幔的黏性对拆沉时间具有更重要影响,岩石圈地幔黏性的减小可促进下地壳榴辉岩的拆沉.若假定下地壳榴辉岩的密度为3.48×10~3kg/m~3(即B_2=1),黏性为5×10~(21)pa s,当岩石圈地幔的黏性减小到2×10~(21)pa s时,60km×30km范围的下地壳榴辉岩的拆沉时间可小于18Myr.     

地幔柱与岩石圈相互作用过程的数值模拟

地幔柱的研究是地球科学研究的热点之一.本文主要集中研究地幔柱与岩石圈的相互作用过程.基于质量守恒方程、动量守恒方程和能量守恒方程,通过有限元数值方法可以计算得到地幔柱与岩石圈相互作用的温度场、速度场和有效黏度等的时空图.本文的流变本构模型主要基于非牛顿流体的有效黏度模型,通过数值模拟计算分析了地幔柱与岩石圈相互作用过程,着重讨论了地壳流变结构对此过程的影响.数值模拟结果显示,地幔柱与岩石圈的相互作用分为三个阶段:地幔柱上升期,时间持续到0.2 Ma,平均速度为2.75m·a-1,地幔柱顶部地形开始向上隆起;地幔柱与岩石圈纵向作用期,时间从0.2 Ma到0.26 Ma,地幔柱上升的平均速度为0.83m·a-1,地表地形隆升达到最大值;地幔柱与岩石圈横向作用期,0.26 Ma以后,岩石圈开始剪切变形,地幔柱水平运动速度为0.47m·a-1,当剪切变形达到一定程度,岩石圈底部开始出现拆沉作用.当下地壳流变强度比较小时,上地壳的流变结构控制着地幔柱顶部地表地形隆起程度,流变强度越大,隆升高度越小;而下地壳的流变结构控制着地幔柱两侧地表地形的下沉幅度,下地壳流变强度越小,下沉幅度越大.最后,讨论了数值模拟对峨眉山大火成岩省地幔柱发展演化的应用.     

数值模拟华北克拉通岩石圈减薄的一种可能机制 ——下地壳榴辉岩重力失稳引起的拆沉

通过数值模拟研究,论证了中生代华北拉通岩石圈受挤压,下地壳玄武岩相变为密度较大的榴辉岩,榴辉岩重力失稳能够引起下地壳的拆沉,造成岩石圈的大规模减薄的可能性.通过对比不同规模的榴辉岩减薄方式发现,当榴辉岩规模较大时可能发生双管道拆沉,而规模小时发生单管道拆沉.计算结果表明当榴辉岩的水平尺度为100 km, 200 km, 300 km时,岩石圈拆沉减薄后厚度分别可达92 km,105 km,136 km.对比中生代华北克拉通岩石圈热侵蚀减薄和拆沉减薄机制,从力学机制上都能造成岩石圈减薄,但它们的结果对应着不同的地表地质特征.热对流减薄,中心区域变热变薄,中心区域是一个伸展区域,发展一个大的变质核杂岩区,岩浆的发展时序是从中心向两侧对称分布,并且向两侧减弱拆沉减薄结果中,中心区域受到挤压加厚,反而变冷.在中心区域的两侧发生两个减薄伸展区域,可能对应两个变质核杂岩区域,岩浆的发展时序是从两侧向中心减弱发展.     

大陆不稳定俯冲的动力学研究

利用地球动力学数值模拟方法分析了大陆汇聚速率、地壳流变强度以及岩石圈热结构等三种主要动力学参数对大陆俯冲碰撞模式的影响和制约,并以此为基础,重点对大陆不稳定俯冲的发育条件、主要类型以及动力学演化过程进行探讨.模拟结果表明,有利于弱化大陆岩石圈流变强度的动力学条件,如低流变强度地壳、岩石圈高温热结构或低汇聚速率等,对大陆不稳定俯冲的形成均起到积极作用.根据俯冲大陆的动力学演化特征可将大陆不稳定俯冲概括为"多阶段断离"型、"持续性流入"型和"大规模拆沉"型等三种类型.这三种类型的发育分别与低汇聚速率、中-高速汇聚的大陆岩石圈高温热结构和低汇聚速率下的中等-高流变强度地壳有关.模拟结果也显示,大陆碰撞区域内岩石部分熔融的分布受地壳流变性质的显著影响.这些模拟结果为认知特定地球动力学背景下不同俯冲条件孕育的大陆俯冲格局具有意义.     

秦岭造山带岩石圈电性结构及其地球动力学意义

秦岭造山带三维电性结构特征的详细研究结果，结合其他多种地球物理和地质资料分析发现，自中新生代以来，尤其是晚近时期华北和扬子两地块向秦岭造山带持续陆内俯冲过程中，由于南秦岭岩石圈向北挤入作用，秦岭造山带的后陆冲断褶带和北秦岭厚皮叠瓦逆冲带，现今处于岩石圈叠置加厚与拆沉作用的初始期；与之相反，南秦岭正在经历拆沉－底侵的物质再循环作用，佛坪一带可能发育新的地幔柱；此外，在造山带北、南深部边界与内部不同岩石圈块体之间还伴随强烈的不同性质的走滑作用和物质侧向传输．最后探讨了秦岭造山带构造拆沉与巨厚岩石圈并存以及走滑构造作用等的地球动力学意义．     

上地幔变黏度小尺度对流的数值研究

基于二维模型，利用有限元方法，研究上地幔-岩石圈系统的变黏度小尺度对流. 考虑该系统的黏度随温度以指数形式变化，数值结果表明，当黏度随温度变化较剧烈时，由于低温高黏度，系统的最上部物质不参与对流，系统发育形成一个类似于岩石圈的静止盖层. 计算表面热流、地形起伏及重力异常与对流格局有较好的相关性，高热流、上升地形对应于对流的上升区，反之低热流、下降地形与对流的下降区对应.     

地幔热导率的选取对动力学数值模拟的影响——以岩石圈张裂过程为例

目前存在有多种地幔热导率模型,不同模型在数值和随温压变化的特征上有明显的差异.为探究不同热导率模型对动力学数值模拟结果的影响,本文对不同模型下的岩石圈张裂过程进行模拟研究,探讨地幔热导率对岩石圈热传输、变形和熔融过程的影响及其作用机理.结果显示,不同热导率模型下,岩石圈的变形和熔融特征表现出明显差异.高热导率模型下,岩石圈破裂较晚,形成陆缘较为宽阔,地壳熔融强烈而地幔熔融较弱;低热导率模型下,岩石圈破裂较早,形成陆缘较为狭窄,地幔熔融强烈而地壳熔融较弱.这种差异源于不同地幔热导率下岩石圈和地幔热状态的变化及相应力学性质的改变.高热导率下,热传导的增温效应显著,岩石圈呈现较热的状态,其强度整体较低,壳幔耦合减弱;而低热导率下,热对流的增温效应显著,岩石圈呈较冷的状态,其强度整体较高,壳幔耦合增强.基于模拟结果,本文认为地幔热导率的选取对动力学模拟的结果有着较为显著的影响,相对于随温压的变化,热导率数值的差异对动力学数值模拟的结果影响更大,尤其是对于地幔熔融过程的影响.     

基于P波三重震相研究青藏高原上地幔速度结构及其动力学意义

青藏高原因其复杂的结构和演化历史,一直都是研究大陆碰撞、构造运动及其动力学的热点区域。本文采用三重震相波形拟合技术,基于中国地震观测台网和大型流动台阵记录到的某地震P波垂向记录,获得了包括拉萨、南羌塘和松潘甘孜地块在内的青藏高原上地幔P波速度结构。结果表明:①拉萨和南羌塘地块下方地幔过渡带存在高速异常,推测是俯冲的印度板片滞留体,过渡带底部的板片残余温度较低,使得660-km相变滞后约3～8km。而松潘甘孜地块下方过渡带同样存在高速异常,可能是欧亚岩石圈发生拆沉进入地幔过渡带所致。这说明印度板块俯冲作用的影响已经到达地幔过渡带,其俯冲前缘位于班公怒江缝合带附近。②从拉萨、南羌塘到松潘甘孜地块,200km之上的地幔岩石圈高速盖层速度由南向北逐渐减小,松潘甘孜地块则出现盖层缺失。推测受小规模地幔对流或者热不稳定性的影响,在南羌塘和松潘甘孜地块,增厚的欧亚岩石圈发生拆沉作用,岩石圈被减薄和弱化,造成羌塘地块上地幔低速和松潘甘孜地块上地幔高速盖层缺失。拆沉的冷的欧亚岩石圈可能部分停留在410-km上方,使得410-km抬升约10km,部分沉入地幔过渡带,表现为松潘甘孜地块地幔过渡带中存在高速异常。低温造成660-km下沉约8km,导致地幔过渡带增厚。      

青藏高原现代地壳运动与活动断裂带关系的模拟实验

本文以GPS观测、大地热流测量、较高精度地形数据、全球板块相对运动的REVEL模型为基础,建立了以青藏高原现代构造活动为主要研究对象的东亚地区构造形变场有限元模型.数值模拟结果显示,青藏高原内部和周边地区走滑断裂带的活动对东亚地区地壳运动速率和方向有较大的影响,特别是对青藏高原物质向东南方向运动有显著影响;不同构造块体岩石圈强度的差异直接影响了川滇菱形地块边界断层错动性质.在考虑青藏高原地形附加重力作用和周边板块汇聚作用对现今大型断裂带运动特征控制作用的同时,岩石圈之下的橄榄岩软流圈至转换带物质对流对岩石圈的拖曳力也是必须考虑的底部边界条件.     

1976年唐山地震震时和震后变形的模拟

本文采用三维粘弹性有限元方法拟合唐山地区１９７６-１９８５年观测到的地震震时和震后的水平与垂直地形变，反演华北板块下方深部物质的流变学性质．模型采用多层弹性覆盖层与线性粘弹性层的有限块体，发展断层面上存在着位错运动，并用正交设计法拟合观测数据．模拟计算表明，华北板块下方软流层粘度为７．１×１０^１８Ｐａ·ｓ；上地幔粘度为２．１×１０^１９Ｐａ·ｓ．     

瑞利数对热对流的影响——在地幔柱中的应用

本文分别在直角坐标系和柱坐标系下,研究瑞利数从104逐渐增大到107对热对流的影响,数值计算结果表明:瑞利数越大,地幔柱越窄,地幔柱上升速度也越快;源自上地幔的地幔柱半径的范围为9到210 km.根据峨眉山内带的半径推算出地幔的黏性系数约为3.8×1021Pa·s,地幔柱平均流动速度为2.5 cm/a.     

四川龙门山断裂带高精度地壳/岩石圈黏度结构及其动力学意义

岩石圈黏度是大陆动力学研究中一个重要参数,但是岩石圈黏度,尤其是横向小尺度(<100 km)黏度结构的确定是一个挑战.本文根据电阻率和黏度与它们控制因素的相似关系,直接把一条跨过青藏高原东缘和四川龙门山断裂带的大地电磁(MT)探测的电阻率剖面转换成黏度结构作为输入,在GPS速度和地表地形数据的约束下,利用地球动力学数值模拟获得了该剖面的二维地壳/岩石圈黏度结构.本文推断的黏度与前人获得的区域尺度的黏度值一致,但揭示出了更多的细节.本文的黏度结构揭示出研究区域内的地壳/岩石圈黏度存在较大的空间变化范围(约5量级),黏度值分布在1.48×10^17~8.44×10^22 Pa·s之间;龙门山断裂带下的黏度存在强烈的小尺度横向变化,其中、下地壳的黏度分别为1.99×10^18~8.21×10^20 Pa·s(平均1.17×10^20 Pa·s)和4.09×10^19~7.08×10^20 Pa·s(平均1.77×10^20 Pa·s).基于该黏度结构的地球动力学模型表明驱动青藏高原中-下地壳物质流动的可能是热-化学浮力,以及上地壳和中-下地壳可能处于解耦状态.本文获得的黏度结构可以为龙门山断裂带地震成因和机制、岩石圈小尺度变形和构造应力状态的深入研究提供重要的帮助.     

流变边界层及其对华北克拉通热/地震岩石圈底界差异的意义

二维热传导/对流数值模型显示,纯传导的固体岩石圈与纯对流的流体软流圈之间存在一过渡层,即流变边界层,其间传导与对流共同作用来传递热量.流变边界层厚度主要由软流圈黏性系数（η）控制,而受固体岩石圈厚度及热状态影响很小.随着η从1×10²¹Pa·s降低至1×10¹⁹ Pa·s,流变边界层也随之减薄,流变边界层的厚度与lg（η）成正比. 流变边界层的存在是造成热岩石圈与地震岩石圈厚度差异的重要因素. 全球典型克拉通岩石圈的对比结果表明,地震岩石圈厚度普遍大于热岩石圈厚度,二者的差异多数在70~90 km,很好地验证了流变边界层的存在. 研究发现二者的差异在华北克拉通自西向东逐渐减小：由西部鄂尔多斯的约80 km减少至渤海湾盆地的约20 km. 反映出华北克拉通岩石圈下部流变边界层厚度自西向东减薄,意味着软流圈黏性系数自西向东逐渐降低.这可能与中生代太平洋俯冲脱水形成的低黏大地幔楔有关,从一侧面印证了太平洋俯冲对华北克拉通破坏的影响.     

克拉通岩石圈对流减薄的数值模拟

采用二维有限元数值模拟的方法研究了岩石圈的对流减薄过程,特别是克拉通岩石圈的对流减薄过程.模型的主要参数包括增厚岩石圈的宽度x、增厚倍数γ、以及与岩石圈组分变化导致的黏性和密度变化密切相关的黏性比(ηc)和浮力数(B)或等效密度变化(Δρtc).数值计算结果显示,地幔对流将逐渐减薄增厚的岩石圈部分,(1)当B=0和ηc=1时,即对一般地幔岩石圈,增厚岩石圈对流减薄的时间可表示为0.0073γ0.70 x0.26.将数值结果应用于地球,意味着增厚到300km的岩石圈,如宽度为300km,对流移除增厚部分回到初始平衡厚度120km大约需要225 Ma;如宽度为1500km,移除增厚部分大约需要342 Ma.(2)当B和ηc较小,克拉通岩石圈对流减薄过程与一般加厚岩石圈的对流减薄过程类似,但减薄时间受克拉通组分浮力和黏性比的影响而显著增长,克拉通岩石圈对流减薄的时间可表示为0.0057ηc0.52Δρ-0.21tcγ0.78ηc-0.36 x0.04.因而,对300km厚的克拉通岩石圈,如克拉通岩石圈的密度比周围地幔的密度低0.4%(即B=0.1),宽度1500km,若克拉通岩石圈黏性因组分影响比普通地幔岩石圈大10倍,其被对流减薄到120km大约需要1.18Ga.(3)当B和ηc增大到一定量时(如B≥0.2且ηc10),克拉通岩石圈被移除的过程将发生变化,由于组分浮力的影响,对流主要不是将克拉通岩石圈带到软流圈地幔中,而主要是将较厚的岩石圈物质向两边推送.在此情况下,克拉通岩石圈能长时间(3Ga)保持稳定.     

三维数值模拟研究俯冲区火山的时空活动性

俯冲板块的深部脱水使得上覆地幔含水, 从而降低含水地幔的熔点, 导致上覆地幔部分熔融。 部分熔融的地幔柱一旦喷发到地表就是俯冲带火山, 也形成新的地壳。 相对于周围的地幔来讲, 具有较小密度和黏度的部分熔融地幔的时空活动性就控制着俯冲带火山的时空分布特征。 本文主要回顾近年来运用三维热力学岩石力学模型数值模拟研究与板片脱水相关的俯冲带火山活动的时空分布特性。 结果表明, 部分熔融地幔的有效黏度和密度是影响俯冲板片之上的三维地幔柱横向分布特征的主要因素。 高黏度的部分熔融地幔(10²⁰~10²¹ Pa·s )易于形成近平行于海沟的、 长波长(70~100 km)的、 薄的波状地幔柱; 低黏度(10¹⁸~10¹⁹ Pa·s )的熔融地幔易于形成平行于海沟的, 短波长(30~50 km)的蘑菇状地幔柱和垂直于海沟的山脊状地幔柱。 当部分熔融地幔和周围地幔的密度相差小于50 kg/m³时, 在俯冲板片之上只能形成长波长低幅度(宽50~100 km, 高10~15 km)的地幔山丘。 岩浆产率随着时间的变化反映了火山活动的生命周期性。 板块俯冲速度会影响地幔柱形成的深度和范围大小。 高效率熔融提取会增加新地壳增长总量。 低的板块俯冲速度和低的熔融提取效率会增加上地壳(花岗岩质)和中地壳(英安岩质)化学成分的比例。 数值模拟结果可以很好地解释如日本东北、 新西兰、 南阿拉斯加俯冲区火山的横向分布特征。     

青藏高原下地壳流动方式的数值模拟研究

青藏高原存在柔性下地壳流动被越来越多的学者接受, 但是关于下地壳流动方式及速度存在争议。 地表运动有GPS等直接测量, 上地幔运动有S波分裂间接反映, 下地壳运动目前没有直接观测手段, 使得开展数值分析非常重要。 本文利用三维球壳黏弹性有限元模型研究了青藏高原下地壳柔性流动方式和流动速度。 本文通过对地表GPS观测资料的拟合与不同数值模型的对比分析, 认为青藏高原柔性下地壳东向流动遇到四川盆地的抵阻, 下地壳物质可能仅在高原东南方向存在物质外溢通道, 而在高原东北方向不存在类似的物质通道; 下地壳的流动速度比地表运动速率每年快几毫米至十几毫米, 对应的黏滞系数为10¹⁸~10¹⁹ Pa·s。     

利用震后黏弹性位错理论研究苏门答腊地震(M_w9.3)的震后重力变化

本文利用2003—2011年的GRACE RL05数据提取了苏门答腊地震(Mw9.3)引起的震后重力变化,发现断层两侧震后重力变化速率存在明显差异,断层下盘总体变化率为0.55μGal/yr,断层上盘为0.16μGal/yr.基于子断层叠加的编程思想,本文将Tanaka的黏弹球体位错理论配套计算程序(简称黏弹位错程序)加以改造,克服了其近场计算精度不足(甚至错误)的缺陷,可用来研究大地震引起的近场震后位移与重力变化.本文利用改造后的黏弹位错程序计算了2004年苏门答腊地震(Mw9.3)产生的同震重力变化,计算结果在空间分布和量级上均与利用弹性位错程序计算获得的结果一致,验证了我们对黏弹位错程序进行改造的正确性.最后,结合GRACE卫星观测数据,本文利用Tanaka的黏弹位错理论研究了苏门答腊地区的地幔黏性因子.结果表明,该地区地幔黏滞性具有显著的横向差异,当发震断层上下两盘的地幔黏滞性系数分别取8×1018 Pa·s和1×1018 Pa·s时,模拟的震后重力变化在总体空间分布和变化趋势上与GRACE卫星观测结果更接近.     

A numerical model of exhumation of the orogenic lower crust in the Bohemian Massif during the Variscan orogeny

We present a numerical model of the main phase (370?C335 Ma) of the Variscan orogeny in the central part of the Bohemian Massif. The crustal deformation in our model is driven by radiogenic heating in the felsic lower crust, the lateral contraction of the Moldanubian domain due to convergence with the Saxothuringian plate (in the early stage of orogeny), and the indentation of the Brunovistulian basement into the weakened orogenic root (in the late stage). Our model explains the main geological events inferred from the geological record in the Moldanubian domain: formation of the orogenic plateau and onset of sedimentation at about 345 Ma, rapid exhumation of the orogenic lower crust at about 340 Ma and subsurface flow of crustal material (?? 335 Ma and later). The results of our modeling suggest that delamination of the lithosphere, often invoked to explain the high temperature metamorphism in the orogenic lower crust of the Bohemian Massif, is not the only physical mechanism which can transfer a sufficient amount of heat to the crust to trigger its overturn.