一种基于遗传算法地震震源破裂过程反演方法.docVIP

一种基于遗传算法的地震震源破裂过程反演方法 　　【摘 要】提出一种可以同时采用近/远场地震波反演地震破裂过程的方法，以便可以更好地利用近场高频GNSS数据及远场地震波形数据。在一个给定的区间内以一定步长搜索破裂传播速度，对于每一个破裂速度，采用遗传算法得到归一化的子断层震源时间函数，通过线性反演快速得到各子断层的滑动分布。采用2013年芦山6.6级地震的源-场模型，进行了数值模拟实验，验证了本方法的可行性与稳定性。 　　【关键词】破裂过程反演；遗传算法；高频GNSS 　　0 引言 　　高频GNSS数据已经被越来越广泛地应用于地震学研究。同时，也有越来越多的连续运行GNSS跟踪站被布置到地震带附近，使得我们可以得到许多同震的近场高频GNSS数据。如果仅仅使用这些数据来获得静态同震形变，就会损失数据中的时域信息。本文提出的方法可以将同震高频GNSS时间序列携带了地震震源破裂过程的信息提取出来，以服务于震源物理研究。 　　本文所采用的方法基于遗传算法，将计算过程并行化，并对计算过程进行了重要的分解，提前计算出一些关键的矩阵，能够快速计算每一个“个体”的适应度，从而大大减少了计算所需时间。 　　1 正演 　　将断层面划分为m个子断层。假定纯剪切破裂起始于震中，并以一个（待求的）恒定速度传播，同时假定每个子断层都有统一的震源时间函数（除了破裂开始时间和幅值），每个子断层的开始破裂时间由且仅由破裂传播速度决定，而各子断层震源时间函数的幅值由各子断层的滑动量决定。这个假定也意味着各子断层在破裂过程中不会改变滑动方向。每个子断层的滑动矢量可分解为两个标量：Xm和Ym，其中Xm为滑动矢量在断层走向上的投影，而Ym为滑动矢量在断层倾向上的投影。 　　采用n个上升时间为Δt的斜坡函数来表示各子断层正则化的震源时间函数，其中第一个斜坡函数的幅值a1为1，接下来各斜坡函数的幅值为a2，a3，…，an。 　　台站xj处的理论位移波形由下式表示： 　　Fk（ti，xj）=∑mmanXmfmk（ti-（n-1）Δt-Tm，xj）+∑mmanXmgmk（ti-（n-1）Δt-Tm，xj）（1） 　　其中fmk和gmk为第m个子断层上滑动角为0°和90°的单位滑动对应格林函数的第k分量，Tm为第m个子断层的起始破裂时间，它由破裂速度决定。 　　2 反演 　　其中F观测波形，e为高斯误差。为使求解过程稳定，对各子断层走滑与倾滑分量引入二阶平滑约束阻尼项其中X与Y为为各子断层滑动量组成的矢量，L为矩阵化的离散拉普拉斯算子。 　　在一个合理区间内以给定的步长搜索破裂传播速度，通过使下式取得最小值来得到对应参数X，Y，以及（a2，…，an，…）的最优估计： 　　其中F0为观测波形向量，D为观测波形的协方差矩阵，？茁为经验阻尼系数，Gn为某破裂速度下的格林函数矩阵。 　　这是一个非线性问题。它可以被线性化，通过迭代求解，但这种方法对迭代初始值的选取非常敏感，很容易终止在局部极值处。 　　固定（a2，…，an），将此问题变成一个典型的二次型问题，其最小值由下式得到： 　　上式为一个关于（a2，…，an）的函数，通过遗传算法来得到它的最小值以及对应的参数（a2，…，an）。 　　在遗传算法中，每个候选的矢量（a2，…，an）就是一个“个体”，一定数量的个体组成一个“种群”，而整个遗传、变异、选择、交叉的过程是以种群而非个体为单位。对每个个体，根据（3）式来计算它的适应度，并以此来确定每个个体被选择到下一代的概率，使得“更优”的个体有更大的概率被选择保存下来。 　　由于矩阵Gn过于庞大，若每次将它代入（3）式计算适应度，会使得计算十分缓慢。因此，在遗传算法开始之前，先计算出和，前者是一个小得多的方阵，后者则是一个向量。它们的大小不会随着观测值的增加而变化，而只与子断层数量有关。这一环节至关重要，它大大减小了遗传算法中的计算量，使得本方法有能力使用大量观测数据进行反演。 　　通过遗传算法找到使得（3）式最小的（a2，a3，…）后，便可据下式求得滑动量的最优估计： 　　3 数值实验 　　制作一个检测板，假定一个单位斜坡函数为各子断层统一的震源时间函数，并设破裂从震中开始以3.0km/s的速度传播，采用2013年芦山6.6级地震的源-场模型，选取了26个GPS台站和35个地震台站作为位移/速度波形的接收场点。在合成的波形中，加入RMS为波形平均振幅5%的高斯误差。以此为伪观测波形，对震源破裂过程进行了反演。 　　图1 搜索不同破裂速度得到的波形拟合不符度 　　数值实验对破裂速度从2.6km/s到3.5km/s以步长0.1km/s进行了搜索。在上图中可以看到，波形拟合不符度在“真实”破裂速度3.0km/s时取得最小值。在搜索3.0km