汶川大地震地震波传播的谱元法数值模拟研究
地震波数值模拟方法研究综述
地震波数值模拟方法研究综述在地学领域，对于许多地球物理问题，人们已经得到了它应遵循的基本方程(常微分方程或偏微分方程)和相应的定解条件，但能用解析方法求得精确解的只是少数方程性质比较简单，且几何形状相当规则的问题。
对于大多数问题，由于方程的非线性性质，或由于求解区域的几何形状比较复杂，则不能得到解析解。
这类问题的解决通常有两种途径。
一是引入简化假设，将方程和几何边界简化为能够处理的情况，从而得到问题在简化状态下的解答。
但这种方法只是在有限的情况下是可行的，过多的简化可能导致很大的误差甚至错误的解答。
因此人们多年来寻找和发展了另一种求解方法——数值模拟方法。
地震数值模拟(SeismicNumericalModeling)是地震勘探和地震学的基础，同时也是地震反演的基础。
所谓地震数值模拟，就是在假定地下介质结构模型和相应的物理参数已知的情况下，模拟研究地震波在地下各种介质中的传播规律，并计算在地面或地下各观测点所观测到的数值地震记录的一种地震模拟方法。
地震波场数值模拟是研究复杂地区地震资料采集、处理和解释的有效辅助手段，这种地震数值模拟方法已经在地震勘探和天然地震领域中得到广泛应用。
地震数值模拟的发展非常迅速，现在已经有各种各样的地震数值模拟方法在地震勘探和地震学中得到广泛而有效的应用。
这些地震波场数值模拟方法可以归纳为三大类，即几何射线法、积分方程法和波动方程法。
波动方程数值模拟方法实质上是求解地震波动方程，因此模拟的地震波场包含了地震波传播的所有信息，但其计算速度相对于几何射线法要慢。
几何射线法也就是射线追踪法，属于几何地震学方法，由于它将地震波波动理论简化为射线理论，主要考虑的是地震波传播的运动学特征，缺少地震波的动力学信息，因此该方法计算速度快。
因为波动方程模拟包含了丰富的波动信息，为研究地震波的传播机理和复杂地层的解释提供了更多的佐证，所以波动方程数值模拟方法一直在地震模拟中占有重要地位。
数值地球物理学中的地震波传播模拟
数值地球物理学中的地震波传播模拟地震波传播模拟是数值地球物理学中的重要研究领域，它对于理解地震的机制、预测地震的危险性以及构建地震防灾体系具有重要意义。
本文将从数值地球物理学的角度，探讨地震波传播模拟的原理、方法和应用。
地震波传播模拟是通过计算机模拟地震波在地下介质中的传播过程，以获取地震波的传播路径、传播速度和传播强度等信息。
在地震波传播模拟中，地下介质被划分为离散的网格单元，每个网格单元的物理性质（如密度、速度等）被赋予数值，通过求解弹性波方程，可以模拟地震波在地下介质中的传播。
地震波传播模拟的方法主要有有限差分法、有限元法和谱元法等。
其中，有限差分法是最常用的方法之一。
它将地震波方程离散化为差分方程，通过迭代计算，逐步推进地震波的传播。
有限差分法具有计算效率高、适用范围广等优点，因此在地震波传播模拟中得到广泛应用。
地震波传播模拟在地震学研究中有着广泛的应用。
首先，地震波传播模拟可以帮助我们理解地震的机制。
通过模拟地震波在地下介质中的传播过程，可以揭示地震波与地下介质的相互作用，进而深入研究地震的发生机制和演化规律。
其次，地震波传播模拟对于预测地震的危险性具有重要意义。
通过模拟地震波在地下介质中的传播，可以预测地震波在不同地点的传播路径和传播速度，进而评估地震对人类和建筑物的影响。
这对于制定地震防灾措施、规划城市发展具有重要的指导意义。
此外，地震波传播模拟还可以用于地震勘探和地震监测。
在地震勘探中，通过模拟地震波在地下介质中的传播，可以预测地下构造和地质条件，为油气勘探和矿产资源勘探提供重要信息。
在地震监测中，通过模拟地震波在地下介质中的传播，可以解释地震观测数据，帮助我们更好地理解地震活动的特征和规律。
然而，地震波传播模拟也面临一些挑战和困难。
首先，地震波传播模拟需要大量的计算资源和时间。
由于地下介质的复杂性和地震波的多尺度特性，地震波传播模拟需要高分辨率的模型和大规模的计算，这对计算机的性能和存储能力提出了很高的要求。
地震波波动方程数值模拟方法(严选优质)
地震波波动方程数值模拟方法地震波波动方程数值模拟方法主要包括克希霍夫积分法、傅里叶变换法、有限元法和有限差分法等。
克希霍夫积分法引入射线追踪过程，本质上是波动方程积分解的一个数值计算，在某种程度上相当于绕射叠加。
该方法计算速度较快，但由于射线追踪中存在着诸如焦散、多重路径等问题，故其一般只能适合于较简单的模型，难以模拟复杂地层的波场信息。
傅里叶变换法是利用空间的全部信息对波场函数进行三角函数插值，能更加精确地模拟地震波的传播规律，同时，利用快速傅里叶变换(FFT)进行计算，还可以提高运算效率，其主要优点是精度高，占用内存小，但缺点是计算速度较慢，对模型的适用性差，尤其是不适应于速度横向变化剧烈的模型.波动方程有限元法的做法是：将变分法用于单元分析，得到单元矩阵，然后将单元矩阵总体求和得到总体矩阵，最后求解总体矩阵得到波动方程的数值解；其主要优点是理论上可适宜于任意地质体形态的模型，保证复杂地层形态模拟的逼真性，达到很高的计算精度，但有限元法的主要问题是占用内存和运算量均较大，不适用于大规模模拟，因此该方法在地震波勘探中尚未得到广泛地应用。
相对于上述几种方法，有限差分法是一种更为快速有效的方法。
虽然其精度比不上有限元法，但因其具有计算速度快，占用内存较小的优点，在地震学界受到广泛的重视与应用。
声波方程的有限差分法数值模拟对于二维速度-深度模型，地下介质中地震波的传播规律可以近似地用声波方程描述：)()(2222222t S zu x u v t u +∂∂+∂∂=∂∂ （4-1） (,)v x z 是介质在点（x , z ）处的纵波速度，u 为描述速度位或者压力的波场，)(t s 为震源函数。
为求式（4-1）的数值解，必须将此式离散化，即用有限差分来逼近导数，用差商代替微商。
为此，先把空间模型网格化（如图4-1所示）。
设x 、z 方向的网格间隔长度为h ∆，t ∆为时间采样步长，则有：z∆,i j1,i j +2,i j+1,i j-h i x ∆= （i 为正整数）h j z ∆= (j 为正整数)t n t =∆ (n 为正整数)k j i u , 表示在(i,j)点，k 时刻的波场值。
汶川震源有限差分数值模拟研究
［ 中图分类号］Ｐ １ ． ３ ５３
［ 文献标识码］Ａ
［ 文章编号］ １７ —１０ （０ ２ ０一 ６ ０ ６ ３ ４ ９ ２ １ ）１ Ｎ０ ０— ５
２０ ０ ８年 ５月 １ 日，在 四川 龙 门 山地 区汶 川县 境 内发 生 了 ８ ０级 地震 ，这是 １ ７ ２ ． ９ ６年 唐 山大 地 震 以 来 国 内最 为严重 的 一次地 震灾 害 ，汶川特 大地 震激 发 的地震 波 传播 特 性为 科 学研 究 提 供 了丰 富 的 资料 。
［ 要］２０ 摘 ０ ８年 ５月 １ ２日在 四 川 汶 川 发 生 ８ ０级 地 震 ，震 中位 于 龙 门 山 断裂 带 上 。 地 震 激 发 的地 震 波 可 ．
以用 来 研 究 震 源 破 裂 过 程 及 地 球 内部 的层 次 结 构 。 对 地 震 台站 接 收 到 的地 震 记 录 的 拟 合 图重 新 进 行 了拟 合 ，得 到 一 个 极 易 用 数 学 式 表 达 的震 源 时 间 函 数 ， 并 结 合 龙 门 山 的 实 际 情 况 ， 划 分 为碎 屑 岩 层 、碳 酸 盐 层 及 基 底 层模 拟 区域 ，并 运 用 交 错 网格 的 有 限 差 分 法 对 声 波 方 程 进 行 数 值 模 拟 。研 究 表 明 ，利 用 ３个 不 同初 始 震 源 都 能 得 到 较 好 的波 场 快 照 图 ， 而初 始 震 源 持 续 时 间 短 的 波 场 快 照 精 度 更 高 ； 不 同初 始 震 源 不
根据 均匀 弹性 各项 同性介 质速 度 、 力 、 应 位移 ３个场 变量 以及场 变 量关 系方 程 ， 考虑 剪应 力 的影 响 不 并令 其为 零 ， 二 阶方程 转变 为一 阶方程 ： 将
旦 一
地震波传播模拟中的数值方法
地震波传播模拟中的数值方法一、引言对地球上发生的自然灾害进行研究和预测一直是人类所探究的课题之一。
其中，地震是一种造成极大灾害的自然现象，它的预测和探测对减轻地震对社会影响，提高人类对灾害的应对能力，具有重要意义。
地震波传播模拟是地震研究领域的重要课题，为了更好地预测地震和应对地震灾害，需要对地震波传播的数值模拟方法进行深入研究。
二、地震波传播数值模拟的方法1. 有限差分法（FDTD）有限差分法，英文全称为Finite Difference Time Domain，是一种常用的求解电磁场和声场传播问题的数值方法。
FDTD方法利用有限差分逼近微分算符，将偏微分方程离散化，然后通过差分方程组求解离散化问题。
FDTD方法的优点是较为简便和直观，对于一些基础场问题可以精确求解，但是FDTD方法在离散化问题域时会导致误差，对于具有复杂形状、边界不规则和含有多个介质的问题，其求解需要繁琐的预处理工作和较为复杂的网格划分，求解过程也较为复杂。
2. 有限元法（FEM）有限元法，英文全称为Finite Element Method，是一种广泛应用于工程和科学计算领域的数值方法。
它是通过将一个复杂的问题域分解成多个小问题域，用简单的数学公式在每个小问题域内求解，通过对这些小问题域的求解累加得到整个问题域的解。
FEM方法的特点是能够对不规则的计算域进行处理，求解过程较为直观和简单，对于多介质、弹性、非线性等问题也有很好的处理能力。
但FEM方法对于较为复杂的问题各向异性和自由面的处理比较困难。
3. 间接边界积分法（BEM）边界积分法，英文全称为Boundary Element Method，是近年来发展起来的一种求解偏微分方程的数值方法。
BEM方法将待求解的域分为界面和域外两部分，通过界面上的边界积分求解内部问题。
BEM方法对于不规则和异形问题的边界条件求解有很好的处理能力，并且具有较高的精度和较低的计算量。
但是对于非线性问题处理不够准确，对纯内部问题的求解效果不如其他方法。
ansys汶川地震的发生对周围断层影响的数值模拟研究
汶川地震的发生对周围断层影响的数值模拟研究李玉江，陈连旺，陆远忠，詹自敏中国地震局地壳应力研究所，100085[ 摘要 ] 基于川西-藏东地区三维粘弹性有限元模型，以GPS和震源机制解作为约束条件，利用汶川地震同震静态位移量结果，在考虑地形差异和粘弹性松弛等因素的影响下，研究汶川地震的发生造成周围主要断层库仑破裂应力变化。
研究结果表明：1）汶川地震的发生除造成震源区同震库仑破裂应力减小外，还造成断裂带两端不同程度的应力增加，尤其是断裂带北东段，这与震后较强的余震分布基本吻合。
断层面上库仑破裂应力的变化与静态位移量的分布存在一定的关系，静态位移量的高值区对应库仑破裂应力减小的区域。
2）汶川地震的发生显著增加了鲜水河断裂北西段、东昆仑断裂、龙日坝断裂、岷江断裂以及虎牙断裂库仑应力水平，将提升这些断层面上发震概率。
3）鲜水河断裂带作为菱形块体北部的一条强震活动频发区域，受汶川地震影响较为显著，地震的发生造成断裂带强震复发间隔缩短约52~104a，是未来值得重点关注的强震危险区。
[ 关键词]川西藏东地区，龙门山断裂带，汶川地震，库仑应力变化，数值模拟Numerical Simulation Research on effect of the surrounding faults induced by the Ms8.0 Wenchuan earthquake，Sichuan，ChinaLI Yu-Jiang，CHEN Lian-Wang，Lu Yuan-Zhong，ZHAN Zi-MinInstitute of Crustal Dynamics，CEA，Beijing 100085，China[ Abstract ] In this paper，with GPS and focal mechanism constraints，utilizing the findings of co-seismic static slip distribution，and taking the terrain variance and post-seismicrelaxation effect into consideration，we employed the three dimensional viscoelasticmodel of western Sichuan and eastern Tibet regions to calculate the Coulomb failure stressaccumulation rate on the major active faults and the dynamic stress change caused by theWenchuan earthquake. Our results indicated that：1) The Wenchuan earthquake not onlycaused the Coulomb stress decrease in the source region，but also the stress increase inthe two-terminal，easpecially the northeastern segment，and it’s comparatively consistentwith the aftershock distrbution.Meanwhile，the high concentration area of the static slipdistrbution was corresponding to the Coulomb stress reduction.2)The Coulomb stresschange caused by Wenchuan earthquake showed significant increase on five major faultsegments，which were northwestern segment of Xianshui fault、Eastern Kunlun fault、Longriba fault、Minjiang fault and Huya fault respectively，and these faults need to bepreferentially concentrated in the future.3) As the northern boundary of theSichuan-Yunnan rhombic block，there were frequently strong earthquakes along theXianshuihe fault，also the effect derived from Wenchuan earthquake was distinct，andthe characteristic earthquake recurrence interval of this fault was delayed 52~104a，so itwas ranked the high seismic hazard region in the future.[ Keyword ] Western Sichuan and Eastern Tibet regions，Longmenshan fault zone，Wenchuan earthquake，Coulomb stress change，Numerical simulation1前言2008年5月12日汶川Ms8.0级地震发生在南北地震带中段、青藏高原东缘的NE向龙门山断裂带上。
汶川Ms8.0地震同震效应的三维非线性数值模拟分析
（ ）但是 ，在龙 门山断裂带东北段 ，压性应力增量＞ ２ 张性应 力增量 ，压性背景构造 应力有所加载 ，导致余 震大多发生在断裂带东北段上 ；在龙 门山断裂带西南段 ，压性应 力增量与张性应力增量大致相当，压 性背景构 造应力虽未卸载 ，但是也未加载，可能是断裂带西南段余震较少 的力学机制之一。
位移场 的一致性很强。此种位移场 空间分布 图像既有利 于余震序列发生在龙 门山断裂带东北段，且有利于余震
序列大多为走滑 型地震事件 。 ２ 汶 川地震发生之后，龙门山断裂带应力状态的变化特征
汶川地震主震发生之后 ，余震序列大多发生在龙门山断裂带 的东北段 ，而西南段的余震很少。为什么会有
型地震 。汶川地震主震引起 的位移场 的模拟结果显示，主震破裂发生之后引起 的位移场的空间分布具有如下特 征 ：龙 门山断裂带西南段 由于受到鲜水河断裂带、安宁河断裂带 的控制作用 ，位移矢量方向较杂乱，位移场一
致性较差；而在断裂 带东北段则产生 了几乎平行于断层走 向的位移矢量场 ，垂直于断层走向的位 移分量很小 ，
裂均采用接触摩擦单元处理，整个 ３Ｄ模型划分为单元 ３ ３ 一 ７ ２ ３个，节点 １ ６ ０ ４ ９个 。利用陈运 泰等给 出的汶 川地震主震破裂过程 的研究成果，我们把主震破裂事件分成 ７个破裂子事件 ， 应数值模拟计算的 ７个载荷步。 对 具体模 拟计算步骤如下：第 １ ：施加 由 Ｇ Ｓ观测数据确定的位移边界载荷 ；第 ２ 步 Ｐ ～８步：分别施加对应 ７个
地震波数值模拟与分析
地震波数值模拟与分析地震波是地震活动中最重要的研究对象之一。
而地震波数值模拟和分析则是地震学领域中的重要研究方向之一。
在地震波数值模拟和分析的过程中，人们可以通过计算机模拟地震波的传播过程，并从中获取有关地震特征及其引起的地表破坏和建筑物结构变形等各种信息。
这对于地震灾害的预防、预测和减轻有着重要的意义。
地震波的数值模拟方法主要有有限差分法、有限元法、边界元法和谱元法等。
其中，有限差分法是目前地震波数值模拟中应用最为广泛的一种方法。
有限差分法在解决非线性、多维度和非静态问题方面表现尤为出色。
其基本思想是将地震波场离散成网格，并利用二阶精度差分公式计算各个时刻在网格点处的地震波场值。
有限差分法的优点在于精度高、计算速度快，同时可以对复杂地质构造及其他复杂条件进行模拟分析。
地震波的数值分析方法主要有PTA和TFI等。
其中，PTA是计算地震波传播中频谱组成的一种方法。
PTA方法基于傅里叶变换，将地震波在频域中进行分析，主要考虑波振幅和频率之间的关系。
通过对地震波的频谱进行分析，可以得出波传播路径、应变速率及层间的速度等信息。
而TFI则是通过时间域内的雷克子波分析地震波的能量分布，从而得出地表加速度和地震破坏信息。
当我们研究地震波数值模拟的同时，还要重视地震波分析的意义。
地震波的分析能够帮助我们对地震发生的原因、机制及它们对地表的影响进行研究。
同时，地震波分析也可以帮助我们评估地震对建筑物和基础设施的破坏。
这项工作通常涉及结构动力学模拟、震害评估、震害预测等研究领域。
此外，通过地震波分析，我们也可以了解地震所带来的生态影响和异常现象（如水波、地陷等）。
在地震波数值模拟和分析过程中，实际数据采集十分必要。
地震数据采集主要分为地震观测和近场强动观测两种方法。
地震观测是通过装置地震仪器等方法获得的数据。
而近场强动观测则是通过现场安装观测设备，获取地震波传播的信息。
同时，人工模拟地震波也是一种可行的方法，但其对于地震波的形态和波速等方面需进行较为精确的估计。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
基础 , 基于 AK135 理论地球模型 , 对汶川大地震激 发地震波传播进行数值模拟研究 , 并将数值模拟结 果与实际观测资料进行了拟合对比分析.
.
随着地震波理论在天然地震和勘探地震中的应 用 , 数值模拟地震波传播已成为研究复杂地下结构 的有效方法 , 可以用来探究不同介质构造条件下地 震波场传播特性 . 数值模拟地震波场可以利用数值 的方法 , 重现地震波产生和传播过程中的关键性信 息 . 地震波属于机械波方程, 其数值模拟方法一直在 地震波研究中占有重要地位, 主要包括有限差分法、 有限元法、伪谱法及谱元法等. 有限差分方法是最早 被用于地震波传播数值模拟研究中的 [10,11],能够成功 的应用到非对称结构模型中 , 有效的实现了对横向 不均匀性介质地球模型的地震波传播数值模拟 [12]. 该方法算法简单 , 计算速度快 , 占用内存小 , 但只适 应于相对简单的地质模型 [13]. 伪谱方法又称配置点 法 , 基 于 快 速 Fourier 变 换 、 Chebyshev 变 换 及 Legendre 变换, 被广泛的用于数值模拟区域地震波传 播问题 [14~16] 和全球结构的地震波传播问题 [15], 通过 求解二维柱坐标弹性波动方程结合 Fourier 伪谱法实 现了基于整个地球模型的地震波传播数值模拟 , 也 能 够有效 的考 虑了地 球介 质的横 向不 均匀性
中国科学 D 辑:地球科学
2009 年 第 39 卷 第 4 期: 393 ~ 402
《中国科学》杂志社
SCIENCE IN CHINA PRESS
汶川大地震地震波传播的谱元法数值模拟研究
严珍珍
①②
, 张怀 *, 杨长春 , 石耀霖
关于汶川地震破裂过程的研究结果(图 2(a))来更真实 的再现汶川地震震源的时空特性.
N ijk = N (ε i ,η j , ξ k ) = li (ε )l j (η )lk (ξ ) .
震源项为一双力偶力源, 其一般可以表示为[32]:
f = − M ⋅∇δ ( X − X s ) S (t ),
摘要
汶川大地震震中位于震害活跃的龙门山断裂带上, 其强度超过 1976 年唐山大地震.
关键词
汶川大地震 谱元法 数值模拟 全球地震波 并行计算
地震激发的地震波可以用来研究震源破裂过程及地球内部的层次结构. 从弹性波理论基础出 发, 采用 AK135 理论地球模型, 考虑地表地形、 地球介质衰减及地球椭率等特性, 利用谱元法 结合高性能并行计算, 分别对人工点源和复合源所激发的汶川大地震地震波的全球传播过程 进行数值模拟重现. 数值模拟结果显示了地震波在地球表面的传播形态. 通过比较发现, 复合 源数值模拟结果较点源数值模拟结果更能体现汶川地震震源破裂过程的时空特性. 另外, 将 两种震源的数值模拟结果分别与实际观测台站记录波形资料进行拟合对比, 进一步定性地认 识到多个点源组成的复合源的数值模拟结果较单个点源的数值模拟结果与实际观测资料的拟 合程度更好, 证明了汶川大地震破裂过程为多阶段破裂组成的一个复合破裂过程.
393
.
中国地震信息网资料显示, 有地震记载以来, 震中附 近 200 km 范围内发生过 8 次 7 级以上地震, 最大的 是 1933 年四川茂汶北迭溪 7.5 级地震, 极震区为四川
严珍珍等: 汶川大地震地震波传播的谱元法数值模拟研究
及对今后四川盆地地震灾害进行了初步的危险性评 估
[8,9]
Runge 现象产生虚假波[21]等偏差. 谱元法可以看作是 广义有限元方法的一种 , 是解偏微分方程的一种有 效的数值方法 , 它结合了有限元法处理边界结构的 灵活性和谱方法的高精度快速收敛等优点成为地震 波数值模拟的重要工具
[22]
. Maday 等
[23]
将拉格朗日插
值函数应用于谱元法中, 并与 Gauss-Lobatto- Legendre(GLL)积分相结合, 从而保证对角质量矩阵的产生, 大大简化了算法实现和存储量 . 在达到相同精度前 提下, 谱元法能够通过采用较稀疏单元划分, 很好的 解决了有限元方法因单元内高阶插值带来的 Runge 现象 , 并且易与大规模高性能并行计算数值模拟相 结合
[2~4]
茂汶北迭溪 , 烈度达到Ⅹ度. 中国地震研究专家及地 质灾害研究专家就汶川大地震进行了深刻地探讨和 研究 , 并详细解释了这次地震的破坏性强于唐山地 震的主要原因. 张勇等[5]利用远场台站垂向波形资料 反演得到主震及部分余震的震源机制解 , 解释了地 震破裂时间过程 ; 王卫民等[6]结合了汶川地震区域的 地质构造背景, 基于有限断层震源模型, 利用实际地 震波形记录及 GPS 等观测资料反演得到了地震震源 破裂过程; 张培震等[7]对地震发震断裂的滑动速率、 复发周期及其构造成因进行了研究 , 认为此次大地 震是低滑动速率、 长复发周期和高破坏强度的巨大地 震, 是一种新的地震类型 . 国内外很多专家也对汶川 地震震源机制等特性进行了理论研究 , 并对 GPS 观 测数据进行了分析研究 . 如用层析成像方法反演龙 门山断裂带地下结构变化 , 地震引起的应力变化以
(1)
(2)
其中, M 为地震矩张量, S(t)为震源时间函数, δ (X−Xs) 为狄垃克函数(Dirac Delta)分布, Xs 为震中位置. 谱元法计算中一般是输入 (2)式中的点源地震矩 张量参数, 但考虑到汶川地震的破裂长度 , 我们尝试 分别以复合源和点源作为震源项进行数值模拟. 我们的复合源模型参数依据张勇等人提出的震 源破裂过程模型 , 认为汶川大地震由三个相继发生 的子事件组成. 三个子事件的发震时间、破裂时间函 数和大小均根据张勇等 [5] 的模型参数 ( 图 2(a) 所示 ) 给出 . 复合源需要的地震矩张量解及震源时间函数
.
该方法计算精度高 , 但不能直接处理起伏自由界面 , 容易受到模型结构限制 . 有限元法适宜于模拟任意 地质构造, 可以任意三角形逼近地层界面 , 保证复杂 地层形态模拟的逼真性 [19], 但是低阶有限元法对高 频和短波长信号的模拟效果不好 , 需要提高网格分 辨率来弥补, 计算量大
[20]
. 采用高阶有限元法也会因
( 图 2(b) 所示 ) 参考张勇等研究结果 . 其中 , t1=0.0 s, t2=20.8 s, t3=49.6 s, hd1=5.2 s, hd2=7.8 s, hd3=5.2 s, h1=1.9, h2=6.4, h3=2.4. 我们希望能够利用张勇等[5]
图 2 地震震源时间函数 [5] (a)及本文数值模拟采用的复合 源时间函数 (b)
395
严珍珍等: 汶川大地震地震波传播的谱元法数值模拟研究
⎡ 1.523 × 1028 ⎢ = ⎢ −0.065 × 1028 ⎢ 28 ⎢ ⎣ 0.327 × 10
⎡ M rr ⎢ M 2 = ⎢ Mθ r ⎢ Mφr ⎣ M rθ M θθ M φθ
2008 年 5 月 12 日 14 时 28 分, 四川省汶川县 (31.49°N, 104.11°E)发生 Ms8.0 级大地震 . 这是中国 1949 年以来破坏性最强、波及范围最大的一次地震, 也是 1976 年唐山大地震以来, 国内最为严重的一次 地震灾害, 地震的强度远远超过了唐山大地震, 造成 了巨大的人员和经济财产损失 . 汶川特大地震激发 的地震波传播特性为科学研究本次地震提供了丰富 的资料. 中国及世界各国都对这次大地震进行了迅速报 道和研究[1]. 震后研究结果显示汶川大地震发生在龙 门山断裂带上 , 造成了四川盆地西北缘沉降以及断 层以西部分山区的抬升 . 该区域地质构造复杂
根据张勇等[5]提供的主震矩张量解, 结合本研究 中所采用的坐标系 , 换算得到复合源数值模拟中的 三个主要子事件的地震矩张量解 ( 单位为 dyne cm), 分别为
⎡ M rr ⎢ M1 = ⎢ M θ r ⎢ Mφr ⎣ M rθ M θθ M φθ M rφ ⎤ ⎥ M θφ ⎥ M φφ ⎥ ⎦
[17,18]
1
谱元法及数值模型
谱元法的基本思想是在传统有限元方法的基础
上, 每一个单元内使用谱方法或者伪谱法 , 选取正交 多项式为基的插值基函数 , 在各个单元上通过采用 高阶插值基函数 , 提高插值精度, 并进而提高数值解 的收敛速度 . 利用谱元法求解偏微分方程的步骤类 似于有限元方法 , 在网格划分后 , 进行单元处理时 , 谱元法是将解近似的表示成正交多项式的有限级数 展开式 , 单元基函数一般是用截断的正交多项式展 开来确定 [27,28]. 谱元法实质上是高阶有限元方法 , 一 般选择阶数为 4~10 的插值多项式表示为基函数[26,29]. 本文中我们利用的谱元法结合拉格朗日插值基函数 及 GLL 积分, 保证插值点与积分点在同一位置, 从 而使得生成的单元质量矩阵为对角分布的特性 , 同 时采用了显式求解方法 , 大大简化了算法实现过程 , 且易于实现并行计算 , 使得数值模拟方法变得快捷 , 计 算 时 间 缩 短 , 计 算 效 率 进 一 步 提 高 . Komatitsch 等 [24,26,29] 给出了波动方程理论及利用谱元法数值模 拟地震波传播的详细理论基础. 我们的数值模型采用 AK135 模型[30,31], 将整个 球体划分为 144 层, 考虑了地球介质的衰减特性、地 球椭率效应及地表地形 . 其中对于介质的衰减特性 处理是基于标准线性体 [24,26], 地球椭率根据 Clairaut 方程而得出[21,32], 地形数据为 ETOPO5[33]. 数值模型 的网格单元是 27 结点六面体单元, 将地心作为一个 单独的六面体单元 , 沿此六面体的六个面分别向地 表延伸建立整个模型, 如图 1 所示. 从地心到地表逐 层划分网格, 并且在 ICB, 670 km 及 Moho 面三处分 别进行了网格加密,在 ICB 边界处网格首次进行加密, 加密后的网格数为加密前的一倍, 670 km 处网格数再 次加倍 , Moho 面处网格数第三次加倍 , 以满足接近 地表区域需要更高分辨率的要求. 图 1 为数值模型. 图 1(a)显示了地心、外地核及 地幔, 整个球体由六个大小相等的块体组成; 图 1(b) 显示了整个数值模型形态以及地球地形和震中位置. 本次研究中, 使用的程序是 Jeroen Tromp 等人提