在地球的演化历史中,间或的小行星撞击对地表地质和生物圈演化存在着巨大的影响.如果一个直径1 km的小行星撞击地球,该撞击将会引发强烈的地震效应,产生巨大的海啸(如果撞击发生在海上).如果小行星的直径增加到10 km,其影响可能对整个生物圈都是灾难性的——6500万年前发生在墨西哥湾尤卡坦半岛的小行星撞击,就被广泛地认为是导致K-Pg边界生物大规模灭亡的原因(Alvarez et al., 1980).如果小行星直径再大一些,撞击还能导致对趾点物质的强烈位移.Schultz和Gault(1975)就曾认为在月球和水星表面的巨大撞击盆地对趾点产生的隆起和塌陷,是由撞击产生的P波和汇聚的面波引起的.这些破裂区域因为强度的降低在功能上将类似于岩浆通道从而诱发一些潜在的火山活动(Williams and Greeley, 1994).如果小行星直径足够大,其对动力学的影响可能也是巨大的——在地球形成的早期(4.2—3.8 Ga),撞击就曾对地幔-地壳体系产生过巨大的影响(Glikson,2001).在某些大胆的假设中,撞击甚至被认为是下地幔D″物质的主要来源(Tolstikhin and Hofmann, 2005).
撞击能否对地幔对流格局和地表热流分布也产生明显的影响?Watters等(2009)曾较早地讨论过这个问题.其在模型中考虑了直径为几百千米量级的小行星对地球和火星产生撞击后全地幔对流的变化.撞击的效果被等效为一个大致以撞击最深点为球心,温度沿半径方向衰减的球形热异常.计算结果表明,这种巨大的撞击对地幔对流的影响也是十分巨大的.
虽然有一些证据表明地球曾经经历过直径为几百千米的小行星的撞击(Glikson,2001),但也有数值模拟的结果显示,这种超大的撞击几乎可以使整个行星表面(甚至内部)面目全非(Nimmo et al., 2008),这使得我们很难用一个(或多个)简单的模型来等效撞击效果.因此为使结果更为真实和可信,不同于Watters等(2009),我们在模型中将不考虑这种超大的撞击.
为此,我们考虑了直径为1000 km的陨石坑(小行星直径约100 km)对地幔对流的影响.之所以这样选择,是因为月球、水星、火星上都有着直径1000 km左右的陨石坑,而超过2000 km的却又极少.同时地质学也有着比较好的模型解释这类陨石坑的结构和演化(Head,2010).一般来说,撞击效果在被破坏之前总会以陨石坑的形态存在.陨石坑是有着明确特征的特殊构造(Grieve,1987),稍微大一点的还有着显著的地球物理学特征(尤指重力场和磁场)(Pilkington,1992).陨石坑形成的初期,撞击点周围的物质会迅速升温,部分会熔化,还有一部分会被抛射.如果陨石坑比较大,撞击之后陨石坑底部还会迅速地回升(French,1998;Grieve and Cintala, 1992;Cintala and Grieve, 1998),显然这些因素都会影响到热对流的初始模型.为此我们将撞击效果等效为两种效果的叠加:(1)撞击区的热异常(温度效果),(2)陨石坑底部的回升作用(速度效果).因为我们初步的计算结果表明1000 km的陨石坑还不足以对全地幔对流产生显著的影响,因此我们只集中 考虑了其对上地幔对流的扰动.与Watters等(2009)集中讨论对流在spreading time的状态不同,本文考虑了对流从扰动到新稳定的全过程,以及可能存在的模式和规律.此外,我们还以Chicxulub陨石坑的热模拟结果(Abramov and Kring, 2007)为初始条件,讨论了小行星直径约10 km的较小规模的撞击(形成直径100~200 km的陨石坑)对上地幔对流的影响. 2 模型及参数设置 2.1 直径为100 km小行星的撞击等效模型
Head在其2010年的文章中对多环盆地这类巨型陨石坑的形态及演化给出了较好的解释,如图 1所示,我们将以这个模型为基础建立直径为100 km小行星的撞击等效模型.
从图中可以看出,初始时刻的陨石坑有着很大体积的位移区域和熔化区域,位移区域像一个半球壳,熔化区域像一个半球.此时我们只要知道瞬态坑和融化区域的直径就可以得到陨石坑的初始异常温度分布了.
由Croft(1984)的结果,假定Dtc,DQ,Dr分别是瞬态陨石坑的直径、复杂陨石坑的最小直径(对于地球,复杂陨石坑的最小直径为4 km)、稳态陨石坑的直径,则有如下的关系:
由此得到瞬态陨石坑的直径在313.5~608.7 km 之间,我们取为500 km.而融化区域的直径并不直 接可知.但其体积却可以由公式:
计算,其中Vm是熔化的体积,c和d是系数(Grieve and Cintala, 1992).如果假定撞击速度为15 km·s-1,则c和d分别为6.21×10-4和3.85.同时,如果假定融化区域是一个半球形,则可得熔化半径大致为190 km.由此就可以得到初始异常温度分布.图 2所示分别是模型LV01和LV02(具体参数见后文表 1,表 2)在撞击之后的温度场,可以看出异常温度区域的面积占据了整个对流模型的很大比重.
除此之外,撞击发生后,陨石坑底部还会经历一个迅速的回升过程,并伴随位移区的部分塌陷(French,1998).因为庞大体积的物质迁移仅在短时间内完成,地质学家很难告诉我们回升过程中陨石坑的具体形态是怎样,因此这里我们将采用热对流方程计算出陨石坑回升后的结果.
已知撞击最深点的回弹高度SU和稳态陨石坑的直径D满足关系(Grieve et al., 1981):
可得陨石坑的最大回弹高度约为120 km.同时,如果假定陨石坑的回升过程与冰后回弹类似的话,那么我们可以由冰后回弹的模型类比到陨石坑底部不同时刻的速度场.首先,已知冰后回弹位移随时间变化的公式为:
其中,ω为任意一点在t时刻的回升高度,ωm是该点的最大回升高度,τr是延迟时间,其满足关系:
这里μ,ρ,g,λ分别是黏度,密度,重力加速度和起伏区域的波长(Turcotte and Schubert, 2002),这里我们选用陨石坑直径的2倍,即2000 km作为波长.假定任何一点的最大回升高度由中心点的回升高度SU和该点与中心点的距离x决定,且满足:
则将(6)式代入(4)式,并对(4)式两边关于时间求导,可以计算出陨石坑底部任意一点的回升速度关于时间和位置的关系为:
以此结果限制热对流的速度场,就可以计算出陨石坑回弹以后的温度场.过程如下:
已知Boussinesq近似和无限Pr and tl数下的热对流方程为
其中,ui,p,T,α,gi,κ分别是速度在i方向的分量,动态压强,温度扰动,热膨胀系数,i方向的重力加速度和热扩散系数(Moresi and Solomatov, 1995).如果在某一时间步,温度场已知,我们都可以由(8)和(9)式计算出新的速度场,然后用(7)式在当前时间步的值修改该速度场,最后将修改后的速度场代入(10)式计算出下一时间步的温度场……以此类推.最终所得温度场即为热演化的初始温度场.如图 3所示为LV01和LV02回弹后的温度场.由于冰后回弹的时间很短(根据黏度的不同从几万到几十万年不等),因此,除了陨石坑底部的回升之外,我们得到的新模型跟之前的变化很小(事实上限于我们程序精度的问题,计算结果略有偏差).我们将以这个结果为初始值,研究其热演化的结果和意义.
在考虑直径为10 km的小行星撞击的等效模型时,我们直接选用了Chicxulub陨石坑作为参考.这是因为,一方面产生该陨石坑的小行星直径刚好在10 km左右;另一方面,该陨石坑保存比较完好,同时有过比较全面的地质学和地球物理学勘探(Hildebr and et al., 1991;Pope et al., 1996;Morgan et al., 2000;Melosh,2001;Vermeesch and Morgan, 2004; Vermeesch and Morgan, 2008),资料丰富.
相对于上地幔对流模型,该陨石坑的直径很小.由于程序精度的限制,所以,虽然Chicxulub在形成之后也有过底部回升的现象(Morgan et al., 2000; Melosh,2001;Vermeesch and Morgan, 2004;Vermeesch and Morgan, 2008),但我们不予考虑.而只是简单将撞击效果等效为叠加在稳态温度场上的一个热异常.这里我们主要参考了Abramov和Kring(2007)的数值模拟结果,假定热异常区域的大小为80 km×20 km,异常温度为1500 K.图 4是模型LVS03(撞击点在右侧,黏度为常黏度,模型参数在后文表 1,表 2中)的初始温度场,其中箭头所指的为异常温度区域,可以看出异常温度区域面积远小于整个对流模型.
本文采用二维有限元程序Citcom(Moresi,1995;Leng and Zhong, 2008)计算热对流演化.黏滞流体为不可压缩,计算区域为直角系,上下为自由滑移边界条件,左右为反射边界条件.假定撞击之前上地幔对流已为稳态,同时将一些重要的参数设置为地球的当前值(表 1).虽然有数值模拟的结果表明上地幔与温度相关的黏度会对对流格局产生影响(叶正仁和王建,2003),但此处我们也只是简单地将上地幔黏度设为常数.对于1000 km陨石坑,考虑了黏度为2×1020Pa·s和1021Pa·s时撞击分别发生在热柱右边及热柱上的情况(LV01,LV02,UV01,UV02,见表 2).对于直径10 km的撞击考虑了黏度为2×1020Pa·s时撞击分别发生在热柱左边、热柱上、热柱右边的情况(LVS01,LVS02,LVS03,见表 2).模型的长宽比为3:1,相应的网格点为385和129. 3 结果及分析
(1)对于直径100 km的陨石撞击,如果上地幔的黏度比较低,则撞击之后的热演化有明显的模式 可以遵循,且彼此之间衔接紧密.如图 5所示为LV01 的演化结果(LV01的初始温度场见图 3a).
从图中可以看出,对流会首先经历一个短暂的调整阶段(约50 Ma),之后进入到一个有着四个对流环的亚稳态(约300 Ma),然后会在经历一个相对较长的再调整阶段(约200 Ma)之后达到新稳定.撞击发生在热柱之上(LV02,初始温度场见图 3b)的情况与前者类似,不同的是,在短暂调整(约30 Ma)后,会有一个时间相对较短的六对流环亚稳态(约80 Ma,图 6),之后再是长时间的四对流环的亚稳态(约200 Ma),然后再经历调整(约350 Ma)达到新稳定.从扰动到新稳定,LV01和LV02的用时分别为550 Ma和650 Ma.可见,巨大的撞击会使上地幔的对流发生强烈的变化,同时撞击点的位置对扰动的形态和持续时间也会产生较大的影响.
图 7给出了LV01和LV02的表面热流随时间的关系,同时还在其上标注了不同对流模式所在的时间段.可以看出撞击发生之时的初次调整阶段,热流的波动都会比较大,且没有明显的规律或周期.在多个对流环的亚稳态阶段(包括四个和六个对流环),热流波动呈现较好的周期性,但振幅仍然比较大.最明显的特征出现在稳定之前的最后一次调整,此时波动的周期性开始减弱,热流值会整体小于正常情况下的热流值.热流值减小则内部的平均温度会升高,这似乎表明,对流会先升高自身的温度,变得相对不稳定,然后再通过放热进入一个更低更稳定的态.总体上来讲,对流在这个阶段十分混乱,没有明显规律和形态.之后就是新的稳定阶段,达到稳定后的热流不会立刻变得跟未撞击的热流一样,但只要演化时间够长二者总会趋于一致.
对于高黏度的情况,撞击之后上地幔的热演化模式也有一定的规律,但相较于低黏度情况要简单很多,总体上可以分为两个阶段:调整和新稳定.图 8是UV01和UV02的表面热流图,可以看出,黏度增加后达到新稳态的用时变短了(UV01约500 Ma,UV02约250 Ma).同时,撞击点的位置也会强烈地影响达到稳定所用的时间——撞击点在热柱上的用时明显低于撞击点在热柱、冷柱之间的用时,这个结果刚好和低黏度的情况相反.
(2)对于直径为10km的小行星撞击,在我们考虑的所有情况(LVS01,LVS02,LVS03)中,无论热对流形态还是表面热流,撞击都没能对其产生明显的影响.图 9是LVS03的表面热流图,可以看出热流只在撞击发生后的极短时间(2—3 Ma)内有一些波动,而这还不排除初始时刻引入的计算误差.LVS01和LVS02的情况完全与之类似.由此可见无论撞击点位置如何,直径10 km的陨石撞击对上地幔对流的影响都微乎其微.
(3)对于100 km直径的小行星撞击,热对流稳定以后,热柱都会发生一定程度的位移:对于低黏度的情况,撞击点在热柱上和热柱、冷柱之间时,热柱分别向撞击点移动了0.1和0.3(无量纲,真实值分 别为67 km和201 km);对于高黏度的情况,撞击点在热柱上和热柱、冷柱之间时,热柱分别向撞击点移 动了0.05和0.15(真实值分别为33.5 km和100.5 km). 可见,当撞击点在热柱上时,热柱向着撞击点的移动较小;当撞击点在热柱、冷柱之间时,热柱向着撞击点移动得较多.而位移量的多少与黏度直接相关,黏度越大移动越小.对于10 km的小行星撞击,无论撞击点位置如何,热对流稳定以后,热柱都会向模型的左边移动,这与未发生撞击的正常热演化的结果一致.因此可知,小的撞击不会引起热柱的移动,而我们看到的移动现象,只是计算本身的结果. 4 结论和讨论
本文分别模拟了10 km和100 km直径的小行星的撞击对地球上地幔对流演化的影响.之所以这样选择,是因为这种规模的撞击可以保留比较完整的地质记录.因此,相比于Watters等(2009)基于温压实验把撞击结果等效为一个热异常,我们这种将撞击结果等效为热异常和速度异常叠加的方式或许更接近真实的物理情况.但这也限制了我们研究更大规模撞击的可能,因为地质学很难比较完整地记录小行星直径几百千米的超大撞击.
我们的计算结果表明,直径10 km的小行星撞击对上地幔对流的影响微乎其微,而直径100 km的小行星撞击则足够对上地幔对流格局产生巨大的影响.对于前者,撞击无法引起对流格局的改变,所看到的热柱移动现象也与其无关(即使没有加入热异常,热柱也会移动,这是计算本身的结果).其表面热流的扰动时间大概是2—3 Ma,这与基于热传导的数值模拟结果(1.5—2.3 Ma)(Abramov and Kring, 2007)接近.可见,对于小撞击,热传导的方式足够用来散去异常的热,因此也就无法对对流产生影响.对于后者,撞击不仅强烈地影响对流格局,同时从扰动初始到结束,对流格局的变化还有着很好的模式(特别是低黏度的情况).这可能主要与初始温度场有关.因为初始温度场大致是在稳态温度场上叠加的一个异常温度场,其本身已经很接近一个亚稳态的温度场,所以才会在短暂调整后就再次进入亚稳态的模式.在亚稳态阶段,由于对流环和热柱的变化,模型顶部物质的运动方向也会发生变化.联系到亚稳态的较长持续时间(100—300 Ma),可以料想,对流在亚稳态阶段会影响或改变部分板块物质的运动从而引起板块的裂解或汇聚.与我们研究对流从扰动到新稳态的全过程不同,Watters等(2009)主要研究了对流在spreading time(大致相当于我们模型中的初次调整阶段)的变化及可能的规律.由于其选取的黏度都比较大,因此没有出现与我们低黏度模型中的“调整、多对流环、调整、稳定”模式类似的情况,但却与我们高黏度模型的“调整、稳定”模式十分相近.另一方面,由于没有考虑撞击点位置的影响,其对流在spreading time的变化规律或许也还不够完善.
同时我们知道,整个扰动持续的时间与黏度直接相关,黏度越大时间越短,这是因为黏度越大对流的阻力也越大.撞击点的位置也会明显地影响扰动的持续时间,但二者的关系却不是很清晰:在我们的计算中,低黏度情况下,撞击点在热柱上扰动持续时间较长;而在高黏度情况下,撞击点在热柱上的扰动持续时间反而越短.同时,无论高黏度还是低黏度,巨大的撞击都会引起热柱的位移.位移量与黏度和撞击点都有一定的关系,黏度越大位移越小,同时撞击点越接近热柱位移也会越小.虽然撞击可能会引发强烈的火山喷发(Williams,1994),但是撞击后的表面热流不断衰减,表明它似乎不是地球上超级火山周期性喷发(Mason et al., 2004)的原因,因此超级火山的周期性喷发可能只是由对流本身的不稳定造成的.
本文在讨论小行星撞击对上地幔对流的影响时,假定了上地幔存在稳定的对流、黏度为常黏度.真实的上地幔参与全球对流并耦合小尺度对流,黏度在岩石圈和软流圈之间也会发生强烈变化,所以我们数值模拟的结果需要进一步使用考虑全地幔对流和复杂黏性结构的更为真实的模型加以验证.同时,在考虑100 km小行星撞击的等效模型时,我们的模型更多地采用了地质学结果.事实上,一些模拟撞击的数值结果也可以被采纳,以用来验证我们得到的结果的普适性;或者直接将模拟的撞击结果作为初始条件,进而用来研究直径在几百公里量级的巨型小行星的撞击对全地幔对流的影响. 致谢 感谢冷伟教授对该文的修改建议.
[1] | Abramov O, Kring D A. 2007. Numerical modeling of impact-induced hydrothermal activity at the Chicxulub crater. Meteoritics & Planetary Science, 42(1): 93-112. |
[2] | Alvarez L W, Alvarez W, Asaro F, et al. 1980. Extraterrestrial cause for the Cretaceous-Tertiary extinction. Science, 208(4448): 1095-1108. |
[3] | Cintala M J, Grieve R A F. 1998. Scaling impact melting and crater dimensions: Implications for the lunar cratering record. Meteoritics & Planetary Science, 33(4): 889-912. |
[4] | Croft S K. 1984. The scaling of complex craters. Journal of Geophysical Research: Solid Earth, 90(S02): C828-C842. |
[5] | French B M. 1998. Traces of Catastrophe: A Handbook of Shock-Metamorphic Effects in Terrestrial Meteorite Impact Structures. Houston: Lunar and Planetary Institute. |
[6] | Glikson A Y. 2001. The astronomical connection of terrestrial evolution: crustal effects of post-3.8 Ga mega-impact clusters and evidence for major 3.2±0.1 Ga bombardment of the Earth-Moon system. Journal of Geodynamics, 32(1-2): 205-229. |
[7] | Grieve R A F, Robertson P B, Dence M R. 1981. Constraints on the formation of ring impact structures, based on terrestrial data.// Proceedings of the Lunar and Planetary Science Conference. New York and Oxford, Pergamon Press, 37-57. |
[8] | Grieve R A F. 1987. Terrestrial impact structures. Annual Review of Earth and Planetary Sciences, 15(1): 245-270. |
[9] | Grieve R A F, Cintala M J. 1992. An analysis of differential impact melt-crater scaling and implications for the terrestrial impact record. Meteoritics, 27(5): 526-538. |
[10] | Head J W. 2010. Transition from complex craters to multi-ringed basins on terrestrial planetary bodies: Scale-dependent role of the expanding melt cavity and progressive interaction with the displaced zone. Geophysical Research Letters, 37(2), doi: 10. 1029/2009GL041790. |
[11] | Hildebrand A R, Penfield G T, Kring D A, et al. 1991. Chicxulub Crater: a possible Cretaceous/Tertiary boundary impact crater on the Yucatán Peninsula, Mexico. Geology, 19(9): 867-871. |
[12] | Leng W, Zhong S. 2008. Viscous heating, adiabatic heating and energetic consistency in compressible mantle convection. Geophysical Journal International, 173(2): 693-702. |
[13] | Mason B G, Pyle D M, Oppenheimer C. 2004. The size and frequency of the largest explosive eruptions on Earth. Bulletin of Volcanology, 66(8): 735-748. |
[14] | Melosh J. 2001. Geophysics: Deep down at Chicxulub. Nature, 414(6866): 861-862. |
[15] | Moresi L N, Solomatov V S. 1995. Numerical investigation of 2D convection with extremely large viscosity variations. Physics of Fluids, 7(9): 2154-2162. |
[16] | Morgan J V, Warner M R, Collins G S, et al. 2000. Peak-ring formation in large impact craters: geophysical constraints from Chicxulub. Earth and Planetary Science Letters, 183(3-4): 347-354. |
[17] | Nimmo F, Hart S D, Korycansky D G, et al. 2008. Implications of an impact origin for the Martian hemispheric dichotomy. Nature, 453(7199): 1220-1223. |
[18] | Pilkington M, Grieve R A F. 1992. The geophysical signature of terrestrial impact craters. Reviews of Geophysics, 30(2): 161-181. |
[19] | Pope K O, Ocampo A C, Kinsland G L, et al. 1996. Surface expression of the Chicxulub crater. Geology, 24(6): 527-530. |
[20] | Schultz P H, Gault D E. 1975. Seismic effects from major basin formations on Moon and Mercury. The Moon, 12(2): 159-177. |
[21] | Tolstikhin I, Hofmann A W. 2005. Early crust on top of the Earth's core. Physics of the Earth and Planetary Interiors, 148(2-4): 109-130. |
[22] | Turcotte D L, Schubert G. 2002. Geodynamics. Cambridge: Cambridge University Press. |
[23] | Vermeesch P M, Morgan J V. 2004. Chicxulub central crater structure: Initial results from physical property measurements and combined velocity and gravity modeling. Meteoritics & Planetary Science, 39(7): 1019-1034. |
[24] | Vermeesch P M, Morgan J V. 2008. Structural uplift beneath the Chicxulub impact structure. Journal of Geophysical Research: Solid Earth, 113(B7), doi: 10. 1029/2007JB005393. |
[25] | Watters W A, Zuber M T, Hager B H. 2009. Thermal perturbations caused by large impacts and consequences for mantle convection. Journal of Geophysical Research: Planets, 114(E2), doi: 10. 1029/2007JE002964. |
[26] | Williams D A, Greeley R. 1994. Assessment of antipodal-impact terrains on Mars. Icarus, 110(2): 196-202. |
[27] | Ye Z R, Wang J. 2003. A numerical research on the small-scale convection with variable viscosity in the upper mantle. Chinese J. Geophys. (in Chinese), 46(3): 335-339. |
[28] | 叶正仁,王建.2003.上地幔变黏度小尺度对流的数值研究. 地球物理学报,46(3):335-339. |