论文部分内容阅读
地壳的运动和变形是地壳介质对板块构造和各种动力过程作用的综合响应,蕴含着丰富的地球内部结构、介质特征和深部过程的信息,长期以来地球科学家通过地质学以及传统大地测量学的观测来确定地壳运动和变形特征,难以得到精确的结果。以全球卫星定位系统(GPS)为代表的空间大地测量技术的发展使得一个区域或局部地区尺度上以高效高时空分辨率获取对地壳运动的观测数据成为可能。进入20世纪90年代以后,GPS形变监测网的覆盖面越来越大,GPS测量技术的精度可以达到mm级,因此,GPS技术目前已成为监测地壳运动和地球动力学现象的主要手段。
地壳的同震和震后形变是地震的能量释放在地壳的动力学响应,由于地球介质的复杂性,只能通过数值方法求得。我们可以通过有限元数值方法模拟地震引发的应变和应力场变化,而有了GPS高精度观测资料的约束,我们可以得到可靠的地震动力学模型。
2001年11月14日在昆仑山口西发生了Ms8.1级地震,2008年5月12日,四川汶川发生Ms8.0级地震。这两次地震都是震深都在15km左右的浅源大震,在地表形成了数百公里的破碎带,地表位移明显,并且有高精度GPS观测资料,为我们研究地壳的同震及震后形变提供了很好的条件。
昆仑山口西地震发生在东昆仑断裂带上。该断裂带位于青藏高原东北缘,南侧是巴颜喀拉地块,北侧是柴达木盆地。电法、人工地震资料显示,巴颜喀拉地块的地壳中部电阻率和剪切波速较低,而且莫霍面的深度可能大于柴达木盆地。汶川地震发生在龙门山断裂带上,该断裂带位于青藏高原东南缘,是青藏高原与四川盆地交界处,对该地区以南的临近区域的层析成像和接收函数分析显示,青藏高原东南缘中上地壳具有低剪切波速构造,而龙门山地区是否也是如此,目前还不得而知。
已有研究者做了对2001年昆仑山口西Ms8.1地震的位移场数值模拟,但是他们的目标是正演应变场和应力场,计算地震前后的应变变化和库伦力变化。本文的研究目标是通过试错法来反演地震区域的弹性常数和粘滞系数在横向和纵向上的不均匀性,从地球动力学的角度印证断裂带两侧地壳的力学性质差异。
本文首先划分有限元计算网格并设定外边界条件。模拟昆仑山地震的有限元网格大小为水平1100km×1100km,厚度设为200km,网格节点数为61×41×16;汶川地震的有限元网格大小为水平1900km×1800km,厚度设为200km,网格节点数为62×81×12,边界条件都为地表自由,侧面零位移,底面垂向上没有运动,水平方向自由。
昆仑山地震的内边界条件,即断层位错参数,参考了谭凯Okada模型反演而得的多段破裂结果,汶川地震的的内边界条件参考了USGS的地震台网数据计算结果。对于两次地震,本文都将模拟计算区域从横向和纵向上划分为6到7个部分,每部分的力学常数各不相同,因此每个同震模型对应不同的参数组合。在正演模拟了几千个同震模型后,将模拟的地表位移场与GPS观测数据比较,并结合先验信息,得到了最优的同震模型参数。昆仑山地震的震后模型的形变模式采用Maxwell线性粘弹性体,其力学性质取决于弹性常数和粘滞系数。我们将最优的同震模型的弹性参数组合作为已知条件,通过模拟计算40个震后形变模型在两年内的地表位移场,并与在的GPS在地震11个月后和17个月后的观测结果比较,找出最佳震后模型,探索了粘滞系数的不均匀性。
其中,网格划分和建模的程序以及后处理程序由自己编写,而主计算程序为自己修改后的有限元形变模拟软件TECTON(开源),使用辟节点(split node)技术将断裂面上位移函数的不连续问题转化为内边界条件。
昆仑山地震的同震形变模拟结果表明,东昆仑断裂带南面的巴颜喀喇地块上地壳的杨氏弹性模量约为60GPa,北面的柴达木盆地的上地壳的杨氏弹性模量约为80GPa,而且巴颜喀喇中地壳的泊松比高于柴达木盆地。从整体看来,模拟的震后形变的模式和观测结果大致一致,断裂面南面的震后位移远大于北面的震后位移,反演的巴颜喀喇地块中地壳的粘滞系数比柴达木盆地中地壳的粘滞系数低一个数量级。这些横向上的物性差异解释了同震形变场关于断裂带的不对称,也和电法、人工地震资料相符,意味着青藏高原东北缘的中下地壳可能存在部分熔融。而且,考虑了断裂带南北莫霍面深度差异的模型优于其他模型。
汶川地震的的最优同震模型表明,龙门山地区中地壳的泊松比高达0.30左右,明显高于上地壳(0.25),也高于位于相同深度的成都盆地地壳的泊松比(0.26),符合青藏高原东缘存在中地壳高温低速带的推测。
当然,数值模拟的最优模型的位移场和GPS观测值也有一定的误差,本文对引起误差的各个原因做了详细的讨论和分析,并在最后总结了研究的局限,展望了下一步的工作。