A New Misfit Function for Multimode Dispersion Curve Inversion of Rayleigh Waves
-
摘要: 反演瑞雷波频散曲线能有效地获取横波速度和地层厚度，传统的多模式瑞雷波频散曲线反演需要正确的模式判别.然而，当地层中含有低速软弱夹层或高速硬夹层等复杂结构时，瑞雷波可能会出现"模式接吻"和"模式跳跃"等现象，这些现象极易造成模式误判，进而导致错误的反演结果；同时，传统的频散曲线反演方法需要进行求根运算，进而导致现有的瑞雷波非线性反演速度慢，运算时间长.鉴于此，对传统的Haskell-Thomson频散曲线正演模拟算法进行了改进，提出了一种新颖有效的目标函数.该目标函数直接利用实测频散曲线与迭代更新模型频散函数表面形状进行最佳拟合，无需将多模式频散数据归于特定的模式，可有效避免多模式瑞雷波频散曲线反演模式误识别；同时，该目标函数不需要求根运算，进而大大加快了非线性反演速度.基于粒子群优化算法，利用实际工作中经常遇到的3种典型理论地质模型和某一高速公路路基实测资料进行了理论模型试算和实例分析，检验了本文提出的瑞雷波多模式频散曲线反演新方法的有效性和实用性.Abstract: Inversion of the Rayleigh wave dispersion curve can effectively obtain the shear wave velocity and stratum thickness. And the classical inversion of multimode Rayleigh wave dispersion curve requires correct mode identification. However, there may be "mode-kissing" and "mode-jumping" phenomena in the Rayleigh waves when the stratum contains complex structures such as low-velocity soft intercalations or high-speed stiff sandwich layers. These phenomena can easily lead to mode misidentification, and lead to wrong inversion results. At the same time, the classical dispersion curve inversion method needs to seek the roots, which leads to that the nonlinear inversion of Rayleigh wave is slow and the computation time is long. In view of these, the classical Haskell-Thomson dispersion curve forward modeling algorithm is improved, and a novel and effective misfit function is proposed. The misfit function is directly used to fit the dispersion function surface shape of the iterative updating model by using the measured dispersion curve. It is not necessary to assign the multimode dispersion data to a specific mode, which can effectively avoid the mode misidentification in the inversion of multimode Rayleigh wave dispersion curve. And the misfit function does not require the seeking root operation, thus greatly accelerating the nonlinear inversion speed. In this paper, based on the particle swarm optimization algorithm, three theoretical geological models and a certain roadbed test data often encountered in practical work are used to calculate the theoretical model and analyze the example. And the validity and practicability of the new method of Rayleigh wave multimode dispersion curve inversion are verified.
-
0. 引言
瑞雷波勘探(李庆春等, 2006)具有无损、高效、经济等特点，越来越受到浅地表地球物理和地质工程界的重视，被视为未来最有希望的技术之一.而与该技术密切相关的瑞雷波是一种由纵波和横波干涉形成，并沿着自由表面传播的波.它具有衰减速度慢、信噪比高、抗干扰能力强和在层状介质中的频散特性等特点，已被广泛用于获取浅地表、岩石圈或地幔的横波(S波)速度结构信息(Xia et al., 1999; 胡家富等, 1999; 黄忠贤等, 2014; 李小勇等, 2014; 潘佳铁等, 2014)、岩溶探测(Xia et al., 2004)以及评估浅层土壤的液化特性(Lin et al., 2004)等.在瑞雷波勘探工作流程中，频散曲线的反演是获取地下横波速度信息的关键步骤之一，它主要分为线性化方法与非线性化全局优化方法.线性化的反演方法，如阻尼最小二乘法(Xia et al., 1999)等，不仅严重依赖初始模型的选取，而且需计算偏导数信息，受雅可比矩阵的求取精度的直接影响.而非线性化全局优化算法，如遗传算法(石耀霖和金文, 1995; Dal Moro et al., 2007)、模拟退火法(崔建文, 2004; Pei et al., 2007)、蒙特卡洛算法(Foti et al., 2009)、人工神经网络(Shirazi et al., 2009; 冯杭建等, 2016)、粒子群全局优化算法(Song et al., 2012; 黄发明等, 2015; 谭琨等, 2015)以及小波变换(Tillmann, 2012)等，不仅不严重依赖初始模型的选取，而且不必计算偏导数信息，已被广泛地用于瑞雷波的频散曲线反演中.本文采用的反演优化算法为其中的粒子群全局优化算法，该算法以其实现容易、精度高、收敛快等优点引起了学术界的重视，并且在解决实际问题中展示了其优越性.
在瑞雷波勘探中，当地层中含有低速软弱夹层或高速硬夹层时，瑞雷波频散曲线会出现“之”字形回折现象(张碧星等, 2000, 2002; 刘雪明等, 2009)，且瑞雷波的各模式之间很可能存在“模式接吻”(Liu and Fan, 2012)和“模式跳跃”现象.如果采用传统的瑞雷波频散曲线拟合反演方法，需要对实测数据进行正确的模式判别，而现今对瑞雷波模式判别的研究尚处于理论研究阶段(鲁来玉等, 2006; 刘雪峰等, 2009).Zhang and Chan(2003)分析了错误的模式判别可能产生的影响；Gao et al.(2014)认为非几何波的存在可能导致分析瑞雷波数据处理过程中的模式误判，对原始数据进行切除等预处理可以有效避免这类模式误判的情况发生，但未对更多的典型模型进行分析，所以其解决方案可能不具有普遍性.总而言之，截至目前没有一个切实可靠的方法能对多种地质模型情况下的“模式接吻”或“模式跳跃”现象进行有效判别.如果仅凭主观意识进行模式判别，极易出现模式误判，从而导致错误的反演结果.因此，需要一种有效的方法来解决这一问题.
本文针对上述问题，对面波频散函数正演算法Haskell-Thomson算法(Haskell, 1953)进行了改进，降低了频散函数的数量级，使得频散曲线位于频散函数表面的极小值点处.基于改进后的频散函数表面的这种特性，笔者提出了一种新颖有效的目标函数.该目标函数直接利用实测频散曲线与迭代更新模型频散函数表面形状进行最佳拟合，无需将多模式频散数据归于特定的模式，可有效避免多模式瑞雷波频散曲线反演模式误识别；同时，该目标函数不需要多模式频散曲线求根运算，进而大大加快了非线性反演速度.
本文利用3种实际工作中经常遇到的理论地质模型(两层速度递增型、四层含低速软弱夹层和六层含高速硬夹层)进行了理论模型试算.首先通过瑞雷面波有限差分正演模拟(周竹生等, 2007)得到各模型的理论地震记录，随后利用相移法得到各模型的频散能量谱，并根据能量的极大值来提取瑞雷波频散曲线，以此来模拟实际的瑞雷波勘探工作中获取实测频散曲线的过程.然后基于粒子群优化算法，分别采用传统方法与本文提出的瑞雷波多模式频散曲线反演新方法对模拟的频散曲线进行反演对比分析，检验了新方法的有效性.最后对某一高速公路路基实测资料进行实例分析，检验了新方法的实用性.
1. Haskell-Thomson算法基本原理及其改进
1.1 Haskell-Thomson算法基本原理
Haskell-Thomson算法对层状介质中瑞雷面波的频散函数的计算奠定了利用瑞雷波信号进行地球物理勘探、人工地震勘探、工程勘察、超声无损检测等的基础.后来有多个算法出现，如δ矩阵法、Schwab-Knopoff算法、广义反射-透射系数方法(何耀锋等, 2006)、快速矢量传递算法(凡友华等, 2002).Haskell-Thomson算法的频散函数方程的简化形式为:
$$D\left( {f,c,\mathit{\boldsymbol{m}}} \right) = {\rm{det}}|\mathit{\boldsymbol{U}}\prime \mathit{\boldsymbol{TV}}| = 0,$$ (1) 式中:m代表模型参数；f表示频率；c表示相速度；U′、V为由表面条件(自由表面或埋藏表面)决定的边界矩阵，对于瑞雷波，U′为2×4的矩阵，V为4×2的矩阵:
$${\mathit{\boldsymbol{U}}^\prime } = \left[ {\begin{array}{*{20}{c}} 0&0&1&0\\ 0&0&0&1 \end{array}} \right];\mathit{\boldsymbol{V}} = \left[ {\begin{array}{*{20}{c}} 1&s\\ r&1\\ {2\mu r}&{\mu r}\\ {\mu t}&{2\mu s} \end{array}} \right],$$ (2) 令矩阵Qi=MiPi，其中
$$\begin{array}{l} {\mathit{\boldsymbol{M}}_i} = \left[ {\begin{array}{*{20}{c}} 1&0&0&0\\ 0&1&0&0\\ 0&0&{{\mu _i}}&0\\ 0&0&0&{{\mu _i}} \end{array}} \right];\\ {\mathit{\boldsymbol{P}}_i} = \left[ {\begin{array}{*{20}{c}} 1&1&{{s_i}}&{ - {s_i}}\\ {{r_i}}&{ - {r_i}}&1&1\\ {2{r_i}}&{ - 2{r_i}}&{{t_i}}&{{t_i}}\\ {{t_i}}&{{t_i}}&{2{s_i}}&{ - 2{s_i}} \end{array}} \right], \end{array}$$ (3) 定义传递矩阵Ti=QiEiQi-1，其中
$${\mathit{\boldsymbol{E}}_i} = \left[ {\begin{array}{*{20}{c}} {{{\rm{e}}^{k{r_i}{h_i}}}}&0&0&0\\ 0&{{{\rm{e}}^{ - k{r_i}{h_i}}}}&0&0\\ 0&0&{{{\rm{e}}^{k{s_i}{h_i}}}}&0\\ 0&0&0&{{{\rm{e}}^{ - k{s_i}{h_i}}}} \end{array}} \right].$$ (4) 下面是计算频散函数D(f, c, m)的Haskell-Thomson算法递归公式:
$$\begin{array}{l} {\mathit{\boldsymbol{X}}_1} = {\mathit{\boldsymbol{U}}^\prime },{\mathit{\boldsymbol{X}}_{i + 1}} = {\mathit{\boldsymbol{X}}_\mathit{i}}{\mathit{\boldsymbol{T}}_\mathit{i}}(i = 1,2, \ldots n);\\ D\left( {f,c,\mathit{\boldsymbol{m}}} \right) = {\rm{det}}|{\mathit{\boldsymbol{X}}_l}\mathit{\boldsymbol{V}}|(l = n + 1), \end{array}$$ (5) 式中:i为层数(i=1, 2, 3, …, n)；μi为第i层的硬度，μi=ρiβi2，其中ρi为第i层的密度，βi表示第i层的S波速度；hi为第i层的厚度；ti=2-c2/βi2；k为波数，k=2πf/c；αi表示第i层的P波速度，当c＜αi(c＜βi)时，ri= $\sqrt {\left( {1 - {c^2}/\alpha _i^2} \right)} $ ，si= $\sqrt {\left( {1 - {c^2}/\beta _i^2} \right)} $ ；当c＞αi(c＞βi)时，ri=i $\sqrt {\left( {{c^2}/\alpha _i^2 - 1} \right)} $ =iri，si=i $\sqrt {\left( {{c^2}/\beta _i^2 - 1} \right)} $ =isi.
1.2 Haskell-Thomson算法的改进
在Haskell-Thomson算法中，由于矩阵Ei(公式(4))中含有指数函数，它会引起频散函数的数量级过大，导致频散曲线与频散函数表面的对应关系很不明显，如图 1a所示，表 1为图中频散曲线对应的理论地质模型参数.为此，本文对矩阵Ei进行了改进，在不改变频散函数正负号的情况下降低其数量级.笔者给出了一个因式Ri:
$${R_i} = |{\rm{exp}}\left( {\frac{1}{2}k({r_i} + {s_i}){h_i}} \right)|.$$ (6) 
图 1 Haskell-Thomson算法频散函数表面图与模型A理论频散曲线的对比a.改进前的Haskell-Thomson算法绘制的频散函数表面图；b.改进后的Haskell-Thomson算法绘制的频散函数表面图.图中白色实线为模型A的理论频散曲线Fig. 1. Comparison of the Haskell-Thomson algorithm dispersion function surface and the synthetic dispersion curves of Model A表 1 模型A:三层含低速软弱夹层地质模型参数Table Supplementary Table Model A: a three-layer model with a soft layer trapped between two stiff layers层序号 VS(m/s) VP(m/s) ρ(g/cm3) h(m) 1 220 437 1.8 6 2 160 285 2.0 3 均匀半空间 400 794 2.1 将矩阵Ei除以Ri由此得到改进后的Ei*，即Ei*=Ei/Ri.改进后的传递矩阵Ti*=QiEi*Qi-1，求解频散函数的递归算法与改进前相似.改进后的Haskell-Thomson算法的频散函数表面形状如图 1b所示.对比图 1a和图 1b可以看出，二者的频散曲线是一样的，说明改进后的Haskell-Thomson算法并没有使理论频散曲线位置发生改变，但是频散函数表面的数量级却得到了很大的压制.而且图 1b中的频散曲线与频散函数表面的对应关系很明显，频散曲线均位于频散函数表面的极小值点处.由此可以推断:如果某个理论模型的频散函数表面的极小值点与实测频散曲线对应情况良好，那么该理论模型就能代表实测数据对应的真实地质模型，反之亦然.根据改进后的频散函数表面的这种特性，提出了一种新的目标函数来描述理论模型的频散函数表面与实测频散曲线的对应情况:
$$S\left( \mathit{\boldsymbol{m}} \right) = {\left\{ {{{\sum\limits_{i = 1}^N {\left[ {{w_i}\left| {D\left( {f_i^{{\rm{obs}}},c_i^{{\rm{obs}}},m} \right)} \right|} \right]} }^l}} \right\}^{1/l}},$$ (7) 式中:wi代表第i实测数据点的权重，对于低精度的数据点可分配较小的权重.D(fiobs, ciobs, m)则代表模型为m时，实测数据点(fiobs, ciobs)处的频散函数值，l为范数，在本文中l=2.
在使用新的目标函数时，各数据点是相对独立的，无需将一系列的数据点归于特定模式，故能避免模式误判，并且无需求解频散函数方程的根，省去了大量的计算，大大提高了反演计算速度.
而在传统的曲线拟合方法中，目标函数的形式为:
$$S\left( \mathit{\boldsymbol{m}} \right) = {\left\{ {{{\sum\limits_{i = 1}^N {\left[ {{w_i}\left| {c_i^{{\rm{obs}}} - g\left( {f_i^{{\rm{obs}}},\mathit{\boldsymbol{m}}} \right)} \right|} \right]} }^l}} \right\}^{1/l}}.$$ (8) 传统的目标函数不仅需要进行正确的模式判别，在当发生模式误判时，采用该目标函数可能会得到错误的反演结果.而且在计算频散函数方程的根g(fiobs, m)时会花费很大的计算量，导致反演速度慢.
2. 粒子群优化算法基本原理
粒子群优化算法(particle swarm optimization, PSO)(Poli et al., 2007)首先随机初始化一群粒子(随机解)，在每次迭代过程中，粒子通过跟踪两个极值来更新自己:一个是粒子本身的所找到的最优解，叫做个体极值点(用pbest表示其位置)；另一个是整个种群目前找到的最优解，称为全局极值点(用gbest表示其位置).粒子i的信息可以用D维的向量来表示，位置表示为xi=(xi1, xi2, …, xid…, xiD)T，速度为vi=(vi1, vi2, …, vid…, viD)T，其速度和位置更新方程为:
$$\begin{array}{l} \upsilon _{id}^{k + 1} = \omega \upsilon _{id}^k + {c_1}{r_1}\left( {pbest_{id}^l - x_{id}^x} \right) + \\ {c_2}{r_2}\left( {gbest_{id}^k - x_{id}^k} \right);\\ x_{id}^{k + 1} = x_{id}^k + \upsilon _{id}^{k + 1}, \end{array}$$ (9) 式中:υidk为粒子i在第k次迭代中的第d维的速度；ω为惯性权重；xidk为粒子i在第k次迭代中第d维的当前位置；pbestid是粒子i在第d维的个体极值点的位置；gbestid是整个种群在第d维的全局极值点的位置；r1、r2是[0, 1]上的随机数；c1、c2是加速系数(或称学习因子)，这两项参数分别用于调整粒子的自身经验与社会经验在其运动中所起的作用，表示将每个粒子推向pbest和gbest位置的统计加速项的权重.
在PSO中，用惯性权重ω来控制前面的速度对当前速度的影响，较大的ω可以加强PSO算法的全局搜索能力，而较小的ω能加强局部搜索能力.基本PSO算法是ω=1的情况，因此在后期缺少局部搜索能力，从而导致解的精度不高.在本文中采用了自组织调整惯性权重ω的策略.该策略使得PSO在前期具有较大的ω，从而使PSO具有较好的全局搜索能力；并在后期具有较小的ω，从而使PSO具有良好的局部搜索能力，提高解的精度.自组织的惯性权重定义如下:
$$w{\rm{ = }}\\ \left\{ \begin{array}{l} - {\left( {\frac{{iter}}{{maxiter}}} \right)^2} + {\omega _{{\rm{start}}}}\left( {1 \le iter \le maxiter/2} \right);\\ {\left( {\frac{{iter}}{{maxiter}} - 1} \right)^2} + {\omega _{{\rm{end}}}}\left( {maxiter/2 ＜ iter \le maxiter} \right), \end{array} \right\}$$ (10) 式中:iter为当前的迭代次数；maxiter为最大迭代次数；ωstart和ωend分别为初始的和终止的惯性权重.在本文中，反演时一般设置ωstart为0.9，ωend为0.4.注意:在本文中传统方法与新方法仅采用的目标函数不同，二者的反演优化算法均为PSO.
3. 理论模型试算
3.1 两层速度递增型地质模型
本文首先以一个典型的两层速度递增型理论地质模型(模型B)为例，其模型参数以及反演时的搜索范围如表 2所示.该模型能代表实际中浅地表高速基岩上覆松散沉积物这类地质模型，主要特点为其第1层与第2层的横波速度值相差较大，这容易造成低频段出现“模式接吻”现象.图 2a为利用瑞雷波有限差分法正演模拟得到的模型B的模拟地震记录，图 2b为通过相移法得到的频散能量图.按照图 2b中的频散能量的极大值提取了模拟实测数据，如图 2b中白点所示.当采用传统的曲线拟合的方法对模拟实测数据进行反演时，需要预先进行模式判别，但由于实测数据点几乎处在同一曲线上，没有明显的模式判别标志，因此，极易将所有的实测数据点看作基阶模式.而当采用新方法时是无需进行模式判别的.随后分别采用这两种方法对模拟实测数据进行反演，其结果如图 3所示.在图 3a中尽管两种方法反演得到的模型的频散曲线均与实测数据拟合较好，且在图 3b中二者的第1层的横波速度和厚度也几乎相同，但是第2层的横波速度却相差很大，其中传统方法反演得到的为2 200 m/s左右，明显过大.产生这一差别的主要原因是在13 Hz处存在“模式接吻”现象，模拟实测数据实际上是由10~13 Hz范围内的第一高阶波数据和13~80 Hz范围内的基阶波数据组成的.故实际上传统方法出现了模式误判，导致了错误的反演结果，而新的方法能“自动识别”数据点所在模式，从而得到了正确的反演结果.图 4为模拟实测数据与频散函数表面的对应情况，可以看出，模拟实测数据均位于频散函数表面的极小值点处，且在13 Hz处存在一个低谷(图中黄色线框)，在该处基阶波对应的极小值与第一高阶波对应的极小值发生“模式接吻”，并且模拟实测数据刚好沿着该低谷由基阶波自动过渡到了第一高阶波，从而达到了一种“自动识别”模式的效果.
表 2 模型B:两层速度递增型地质模型参数及反演搜索范围Table Supplementary Table Model B: a two-layer model characterized by S-wave velocities increasing with depth and search space in the inversion层序号 模型参数 搜索范围 VS(m/s) VP(m/s) ρ(g/cm3) h(m) VS(m/s) h(m) 1 150 298 1.8 5 100~300 1~10 均匀半空间 450 802 2.1 ∞ 200~3 000 ∞ 
图 2 基于模型B模拟的理论瑞雷波地震记录及其高分辨率频散能量谱a.60道高精度理论瑞雷波地震记录；b.由图 2a提取的高分辨率频散能量谱.图中白色实点为根据频散能量极大值提取的瑞雷波频散曲线，该频散数据将作为实测频散曲线进行反演Fig. 2. The theoretical Rayleigh wave seismic record and the high-resolution dispersion energy spectra based on model B
图 3 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演模型剖面与模型B的对比(b)图a中:黑色实点为从图 2b中提取的瑞雷波频散曲线；蓝色虚线为利用传统方法反演模型正演模拟的基阶波频散曲线；红色虚点线为利用新方法反演模型正演模拟的多模式频散曲线，该多模式频散曲线在13 Hz处相速度几乎一样，此现象称为“模式接吻”现象.图b中:绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型B；蓝色虚线表示传统方法反演获得的横波速度剖面；红色虚点线表示新方法反演获得的模型剖面Fig. 3. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a), and comparison of the S-wave velocity profiles inverted by classical method and new method and the profile of model B (b)
图 4 新方法反演模型的频散函数表面与实测频散曲线的拟合对比Fig. 4. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve另外，在进行反演工作之前，模型的搜索范围的选择是需要慎重考虑的.不适当的模型搜索范围很可能会使反演算法陷入局部极小值，无法收敛到全局极小值，从而得到错误的反演结果.为此，对上述两种方法对应的目标函数表面形状进行了分析.将厚度设为定值，第1层与第2层的横波速度分别设为横纵坐标，来分析目标函数随着第1层与第2层速度值的变化，如图 5所示.图 5a和5b分别为传统方法和新方法对应的目标函数表面形状，从图 5a中可以看出，在第1层速度为150 m/s左右的位置存在一个狭长的低谷，它能较好地约束第1层的速度，但对第2层的速度的约束能力不强，而且当半空间横波速度增加时，其目标函数值也确实是在减小，使得传统的反演方法得不到正确的结果.而在图 5b中新方法对应的目标函数表面要比传统的复杂得多，它存在多个极小值，正确的解对应其右边部分的全局极小值，还有其他非正确解对应着其左边部分的局部极小值.由于局部极小值的存在，在使用新方法的时候就需要谨慎选择模型搜索范围:一般将第1层的速度搜索下限值设置为靠近第1层真实的横波速度，其中真实的横波速度可以通过观察高频频散能量对应的相速度估测出来，从而使得模型搜索范围尽量避开图 5b中左边部分的局部极小值，以此来避免反演算法陷入局部极小值，从而有利于模型参数快速收敛到全局最优解.在下文的几个例子中，笔者也采用了类似的设置模型搜索范围的策略.
图 5 两种方法对应的目标函数表面形状a.传统方法的目标函数表面形状；b.新方法的目标函数表面形状.二者中的红色圆点代表真实模型，绿色圆点代表各自反演得到的模型Fig. 5. The misfit function surfaces of two methods在此例中，统计了一下两种方法各自花费的时间，在数据量和反演参数完全相同的情况下，采用传统方法所花时间是新方法的102.8倍.由此可见新的反演方法不仅能避免模式误判，还能在很大程度上提高反演计算速度.
3.2 四层含低速软弱夹层地质模型
接下来本文以一个含低速软弱夹层的四层理论地质模型(模型C)为例，其模型参数以及反演时的搜索范围如表 3所示，图 6a为该模型的理论地震记录，图 6b为该地震记录经相移法变换得到的频散能量谱，按照图 6b中的频散能量的极大值提取了模拟实测数据，如图中白点所示.当采用传统的曲线拟合的方法对模拟实测数据进行反演时，由于实测数据点几乎处在同一曲线上，没有明显的模式判别标志，极易将所有的实测数据点均看作基阶模式，而当采用新方法时是无需进行模式判别的.笔者分别采用这两种方法对该模拟实测数据进行反演，其结果如图 7a和图 7b所示.在图 7a中两种方法反演得到的模型对应的频散曲线均与实测数据拟合的较好.在图 7b中新方法反演的模型与真实模型比较符合，其对应的频散函数表面也与模拟实测数据对应很好，模拟实测数据正好位于频散函数表面的最低低谷处，如图 8所示.但是传统方法反演的模型与真实模型差别较大，主要体现在传统方法反演的第2层的速度相对较低，厚度相对较大，这主要是由于在30 Hz和58 Hz处存在“模式接吻”现象，导致传统方法出现了模式误判，产生了错误的反演结果.另外，传统方法反演的半空间的速度也出现了较大的偏离，其原因可能是:在传统方法中，低频数据对半空间的横波速度敏感性相对较高，而且从频散能量谱上拾取的实测数据难免存在误差.低频数据，尤其是最低低频点的较大拾取误差会引起传统的目标函数值发生较大变化，从而引起半空间速度的较大偏离；而采用新的方法时，最低低频点处在频散函数表面的低谷处，当该点发生上下偏移时，其仍处在低谷附近，不会引起新的目标函数的较大变化，相应地也不会引起半空间横波速度发生较大的偏离.由此可见，在这种情况下，采用新的目标函数还可以避免因低频数据的拾取误差造成的对半空间横波速度的过高估计.
表 3 模型C:四层含低速软夹层地质模型参数及反演搜索范围Table Supplementary Table Model C: a four-layer model with a soft layer trapped between two stiff layers and search space in the inversion层序号 模型参数 搜索范围 VS(m/s) VP(m/s) ρ(g/cm3) h(m) VS(m/s) h(m) 1 200 490 2.0 5 120~300 2~8 2 160 392 2.0 5 100~300 2~8 3 260 637 2.0 5 100~500 2~8 均匀半空间 380 931 2.0 ∞ 200~800 ∞ 
图 6 基于模型C模拟的理论瑞雷波地震记录及其高分辨率频散能量谱a.60道高精度理论瑞雷波地震记录；b.由图 6a提取的高分辨率频散能量谱.图中白色实点为根据频散能量极大值提取的瑞雷波频散曲线, 该频散数据将作为实测频散曲线进行反演Fig. 6. The theoretical Rayleigh wave seismic record and the high-resolution dispersion energy spectra based on the model C
图 7 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演模型剖面与模型C的对比(b)图a中:黑色实点为从图 6b中提取的瑞雷波频散曲线；蓝色虚线为利用传统方法反演模型正演模拟的基阶波频散曲线；红色虚点线为利用新方法反演模型正演模拟的多模式频散曲线.图b中:绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型C；蓝色虚线表示传统方法反演得到的模型剖面；红色虚点线表示新方法反演获得的模型剖面Fig. 7. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a) and comparison of the S-wave velocity profiles inverted by classical method and new method and the profile of model C (b)
图 8 新方法反演模型的频散函数表面与实测频散曲线的拟合对比Fig. 8. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve在此例中，两种方法反演的数据量与反演参数完全相同的情况下，采用传统方法所花时间是新方法的149.2倍.进一步可见新的反演方法能在很大程度上提高反演计算速度，而且对于不同的模型与数据，其计算速度的提高程度也是有一定差别的.
3.3 六层含高速硬夹层地质模型
最后本文以一个含高速硬夹层的六层理论地质模型(模型D)为例，其模型参数以及反演时的搜索范围如表 4所示，图 9a为该模型的理论地震记录，图 9b为该地震记录经相移法变换得到的频散能量谱，图中明显存在两个“模式跳跃”现象，分别发生在26 Hz和38 Hz处，频散能量从一个模式跳转至另外一个模式.图中频散能量大致分为3个部分，将其编号.按照图 9b中的频散能量极大值提取模拟实测数据，如图中白点所示.当采用传统的曲线拟合的方法对模拟实测数据进行反演时，根据能量团1、3的走势，极易将二者整个视为基阶波，并同时将能量团2视为第一高阶波.仍分别采用两种方法对该模拟实测数据进行反演，其结果如图 10a和图 10b所示.在图 10a中，两种方法反演得到的模型的频散曲线均与实测数据拟合得较好.在图 10b中新方法反演的模型与真实模型比较符合，其对应的频散函数表面也与模拟实测数据对应很好，如图 11所示；但是传统方法反演的模型与真实模型差别较大，主要体现在传统方法反演的模型的第2层的速度过低，厚度过薄，而且第3层的速度过大，厚度也较大.这主要是由于传统方法出现了模式误判，将实际为第一高阶的能量团3视为基阶波，并同时将实际为第二高阶的能量团2视为第一高阶波，从而得到了一个实质性的偏离了真实模型的错误地质模型.
表 4 模型D:六层含高速硬夹层地质模型参数及反演搜索范围Table Supplementary Table Model D: a six-layer model with a stiff layer sandwiched between two soft layers and search space in the inversion层序号 模型参数 搜索范围 VS(m/s) VP(m/s) ρ(g/cm3) h(m) VS(m/s) h(m) 1 200 490 2.0 2.5 150~400 0.1~5 2 260 637 2.0 2.5 10~400 0.1~5 3 120 294 2.0 2.5 50~600 0.1~5 4 240 588 2.0 3 100~600 0.1~5 5 260 637 2.0 3 150~600 0.1~5 6 370 906 2.0 ∞ 200~600 ∞ 
图 9 基于模型D模拟的理论瑞雷波地震记录及其高分辨率频散能量谱a.60道高精度理论瑞雷波地震记录；b.由图 9a提取的高分辨率频散能量谱.图中白色实点为根据频散能量极大值提取的瑞雷波频散曲线，该频散数据将作为实测频散曲线进行反演Fig. 9. The theoretical Rayleigh wave seismic record and the high-resolution dispersion energy spectra based on the Model D
图 10 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演模型剖面与模型D的对比(b)图a中:黑色实点为从图 9b中提取的瑞雷波频散曲线；蓝色虚线为利用传统方法反演模型正演模拟的多模式瑞雷波频散曲线；红色虚点线为利用新方法反演模型正演模拟的多模式频散曲线.图b中:绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型D；蓝色虚线表示传统方法反演得到的模型剖面；红色虚点线表示新方法反演获得的模型剖面Fig. 10. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a), and comparison of the S-wave velocity profiles inverted by classical method and new method and the profile of model D (b)
图 11 新方法反演模型的频散函数表面与实测频散曲线的拟合对比Fig. 11. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve在此例中，由于模式误判，传统方法反演得到的模型与真实的模型的各个参数的平均相对误差高达77%.而新方法能避免模式误判，其平均相对误差仅为4%，相对于传统方法提高了73%的反演解释精度.
由此可见，该方法能避免模式误判，从而大大提高反演解释精度.
3.4 随机噪声以及层数和厚度的设置对新方法反演结果的影响
在实际应用中，由于实测地震记录不可避免地会含有噪声，从而引起数值计算结果误差.因此，检验新方法的抗噪能力是有必要的.以文中的模型B为例，由于瑞雷波的能量较强，在高精度理论瑞雷波地震记录(图 2a)中分别加入了噪声水平较高的20%的随机噪声(图 12a)，以及40%的随机噪声(图 12c).图 12b和图 12d分别是图 12a和图 12c采用相移法提取得到的高分辨率能量谱，由图可见，加入了20%随机噪声的能量谱与不含随机噪声的能量谱(图 2b)相比，其有效频带范围明显有所减小，尤其体现在高频段部分，有效频带范围从无噪声的10~80 Hz缩小到了含20%随机噪声的10~44 Hz，加入了40%随机噪声的能量谱的有效频带范围也进一步缩小到了11~36 Hz.根据二者的有效频带范围来提取频散曲线，如图 12b和12d中的白色实心点所示，然后采用新的方法对其进行反演，反演结果如图 13所示.不含噪声的数据反演结果与真实模型是最吻合的，其平均相对误差仅为4.4%，含20%的随机噪声的数据反演结果次之，其平均相对误差为5.1%，而含40%的随机噪声的数据的反演模型剖面与真实模型剖面出现了一定的偏离，其平均相对误差为7.7%.由此可见，当噪声水平增大时，其反演结果的误差也随之增大，但三者的反演模型剖面均未实质性地偏离真实的模型剖面，说明新方法对随机噪声具有一定的容忍度，具有良好的抗噪能力.
图 12 基于模型B模拟的含随机噪声的瑞雷波地震记录及其高分辨率频散能量谱a.含20%随机噪声的地震记录；b.由图a提取的高分辨率能量谱；c.含40%随机噪声的地震记录；d.由图c提取的高分辨率能量谱；图b和d中白色实心点为根据频散能量极大值提取的瑞雷波频散曲线，该频散数据将作为实测频散曲线进行反演Fig. 12. Rayleigh wave seismic records with random noise and their high-resolution dispersion energy spectra based on model B
图 13 含不同随机噪声数据反演模型剖面与模型B的对比图中绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型B；红色实线表示无噪声数据反演得到的模型剖面；洋红色点划线表示含20%随机噪声数据反演获得的模型剖面.蓝色点虚线表示含40%随机噪声数据反演获得的模型剖面Fig. 13. Comparison of the S-wave velocity profiles inverted by different random noise data and the profile of model B在本节理论模型试算中，都是已知层数来进行反演的，而在实际的实测资料反演时，层数可能并不总是已知的.因此，对反演时不同层数设置对计算结果的影响进行了探讨.以模型C对应的模拟实测数据为例，当层数设置为4层时，其反演结果的平均相对误差为5.1%.为了提高反演模型的分辨率，将层数设置为5层，反演结果如图 14中洋红色点划线所示，反演模型的平均相对误差为8.5%.为进一步提高反演模型的分辨率，将层数设置为6层，其反演结果如图 14中蓝色点虚线所示，反演模型的平均相对误差为11.9%.由此可见，随着层数的增加，层厚的减小，反演模型的分辨率随之提高，反演模型的误差也会随之增大，但三者的反演模型均未实质性地偏离真实的模型，说明新方法对于不同的层数剖分和厚度设置均具有良好的适用性.另外，这些模拟结果也表明，并不是将地下介质模型剖分的层数越多越好.在实际地球物理资料反演时，为了获得稳定的反演结果，必须在模型分辨率和反演误差之间进行折中选取.
图 14 不同层数设置时反演模型剖面与模型C的对比图中绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型C；红色实线表示设置层数为四层时反演得到的模型剖面；洋红色点划线表示设置层数为五层时反演获得的模型剖面.蓝色点虚线表示设置层数为六层时反演获得的模型剖面Fig. 14. Comparison of the S-wave velocity profiles with different layer setting and the profile of model C4. 实例分析
该实例来自中国河南一条高速公路的路基，路基表面已进行了振动碾压，且表面地形是平的，没有明显的地形起伏.为了查明路基中的低速软弱地层，进行了瑞雷波检测实验.选用的采集参数为:观测方式采用多次覆盖观测系统，仪器采用国产SE2404数字地震仪；4.5 Hz垂直分量低频检波器接收，锤击激发；记录道数为24道，记录长度为0.5 s，采样间隔为0.5 ms, 采样点数为1 024；偏移距为7.0 m, 道间距为1.0 m.图 15a为采集到的地震记录，图 15b为将该地震记录进行f-k变换后得到的f-k域频散能量谱.在f-k域中根据能量极大值提取频散曲线，并将频散曲线变换到f-v域，得到需要进行反演的瑞雷波频散曲线(图 16a中黑色圆点).分别采用本文中的两种反演方法对该实测相速度曲线进行反演，将其分为5层，然后根据实测数据点的分布规律来划定反演时各个地层参数的搜索范围:设置第1层横波速度搜索范围为150~400 m/s，其他层横波速度搜索范围均为100~400 m/s，各层厚度参数搜索范围均为1~8 m.当采用传统方法时，由于实测数据没有明显的模式判别标志，极易将整个实测数据视为基阶波数据来进行反演.两种方法的反演结果如图 16和图 17所示，在图 16a中，两种方法的反演模型的频散曲线均与实测数据拟合较好，在图 17中，其对应的频散函数表面也与实测数据对应得很好.但是将二者得到的横波结构剖面与钻孔资料进行比对后，如图 16b所示，可以发现，新方法反演得到的模型剖面与钻孔资料吻合较好，能明显地揭示出地层中危险的低速软弱夹层，而在传统方法反演的模型剖面中几乎不能揭示这一低速软弱夹层的存在.产生这种差别的原因是实测数据中实际含有多个模式，而非单一的基阶模式，存在着“模式接吻”现象，如图 16a中绿色线框所示.这个现象导致了模式误判，使得传统的反演方法得到了错误的结果.
图 15 某高速公路路基瑞雷波勘探实例a.野外实测的24道瑞雷波地震记录；b.由图 15a中的多道瑞雷波炮集地震记录提取的f-k域频散能量谱.图中白色实线为提取的f-k域实测瑞雷波频散曲线Fig. 15. Exploration case of Rayleigh wave of a highway roadbed
图 16 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演的横波速度剖面与钻孔资料的对比(b)图a中:黑色实点为从图 15中提取得到的实测频散数据变换到f-v域的频散曲线; 红色虚点线为利用新方法反演模型正演模拟获得的多模式频散曲线；蓝色虚线表示利用传统方法反演模型正演模拟获得的基阶波频散曲线.图b中:绿色虚线为粒子群优化算法反演时模型参数搜索范围；带有圆圈的实线表示钻孔资料剖面；红色虚点线表示新方法反演的横波速度剖面；蓝色虚线表示传统方法反演得到的横波速度剖面Fig. 16. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a), and comparison of the S-wave velocity profiles inverted by classical method and new method and the borehole data (b)
图 17 新方法反演模型频散函数表面与实测频散曲线对比拟合情况Fig. 17. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve5. 结论
传统的多模式频散曲线反演需要将多模式频散数据分别归于特定的模式，因而极易出现多模式误识别，导致错误的地质解释.同时，在正演模拟过程中需要求根运算，导致传统的多模式频散曲线反演速度慢，计算时间长.本文针对上述不足，对传统的Haskell-Thomson算法进行了改进，提出了一种新颖有效的目标函数.理论模型试算与实测资料分析结果表明:
(1) 本文利用提出的新的目标函数进行多模式瑞雷波频散曲线反演时，无需将多模式频散数据点归于某一特定模式，即无需进行模式判别，从而有效地避免了多模式频散曲线反演极易出现的模式误判这一技术难题，大大提高了多模式瑞雷波反演解释精度，以文中模型D为例，其反演解释精度提高了73%.
(2) 本文提出的新的目标函数相对于传统目标函数进行多模式瑞雷波频散曲线反演时无需进行频散曲线求根运算，从而大大提高了多模式瑞雷波非线性全局优化反演的计算速度，有效节省了计算时间，以文中模型C为例，其计算速度提高了149.2倍.
(3) 本文提出的多模式瑞雷波频散曲线反演新方法不仅能避免模式误判，提高反演解释精度，而且计算速度快，可有效地应用于多模式瑞雷波频散曲线非线性研究.
致谢: 感谢两位匿名审稿专家对本文进行了细致的审阅并提出了宝贵的修改意见，在此表示衷心感谢！ -
图 1 Haskell-Thomson算法频散函数表面图与模型A理论频散曲线的对比
a.改进前的Haskell-Thomson算法绘制的频散函数表面图；b.改进后的Haskell-Thomson算法绘制的频散函数表面图.图中白色实线为模型A的理论频散曲线
Fig. 1. Comparison of the Haskell-Thomson algorithm dispersion function surface and the synthetic dispersion curves of Model A
图 2 基于模型B模拟的理论瑞雷波地震记录及其高分辨率频散能量谱
a.60道高精度理论瑞雷波地震记录；b.由图 2a提取的高分辨率频散能量谱.图中白色实点为根据频散能量极大值提取的瑞雷波频散曲线，该频散数据将作为实测频散曲线进行反演
Fig. 2. The theoretical Rayleigh wave seismic record and the high-resolution dispersion energy spectra based on model B
图 3 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演模型剖面与模型B的对比(b)
图a中:黑色实点为从图 2b中提取的瑞雷波频散曲线；蓝色虚线为利用传统方法反演模型正演模拟的基阶波频散曲线；红色虚点线为利用新方法反演模型正演模拟的多模式频散曲线，该多模式频散曲线在13 Hz处相速度几乎一样，此现象称为“模式接吻”现象.图b中:绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型B；蓝色虚线表示传统方法反演获得的横波速度剖面；红色虚点线表示新方法反演获得的模型剖面
Fig. 3. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a), and comparison of the S-wave velocity profiles inverted by classical method and new method and the profile of model B (b)
图 4 新方法反演模型的频散函数表面与实测频散曲线的拟合对比
Fig. 4. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve
图 5 两种方法对应的目标函数表面形状
a.传统方法的目标函数表面形状；b.新方法的目标函数表面形状.二者中的红色圆点代表真实模型，绿色圆点代表各自反演得到的模型
Fig. 5. The misfit function surfaces of two methods
图 6 基于模型C模拟的理论瑞雷波地震记录及其高分辨率频散能量谱
a.60道高精度理论瑞雷波地震记录；b.由图 6a提取的高分辨率频散能量谱.图中白色实点为根据频散能量极大值提取的瑞雷波频散曲线, 该频散数据将作为实测频散曲线进行反演
Fig. 6. The theoretical Rayleigh wave seismic record and the high-resolution dispersion energy spectra based on the model C
图 7 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演模型剖面与模型C的对比(b)
图a中:黑色实点为从图 6b中提取的瑞雷波频散曲线；蓝色虚线为利用传统方法反演模型正演模拟的基阶波频散曲线；红色虚点线为利用新方法反演模型正演模拟的多模式频散曲线.图b中:绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型C；蓝色虚线表示传统方法反演得到的模型剖面；红色虚点线表示新方法反演获得的模型剖面
Fig. 7. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a) and comparison of the S-wave velocity profiles inverted by classical method and new method and the profile of model C (b)
图 8 新方法反演模型的频散函数表面与实测频散曲线的拟合对比
Fig. 8. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve
图 9 基于模型D模拟的理论瑞雷波地震记录及其高分辨率频散能量谱
a.60道高精度理论瑞雷波地震记录；b.由图 9a提取的高分辨率频散能量谱.图中白色实点为根据频散能量极大值提取的瑞雷波频散曲线，该频散数据将作为实测频散曲线进行反演
Fig. 9. The theoretical Rayleigh wave seismic record and the high-resolution dispersion energy spectra based on the Model D
图 10 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演模型剖面与模型D的对比(b)
图a中:黑色实点为从图 9b中提取的瑞雷波频散曲线；蓝色虚线为利用传统方法反演模型正演模拟的多模式瑞雷波频散曲线；红色虚点线为利用新方法反演模型正演模拟的多模式频散曲线.图b中:绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型D；蓝色虚线表示传统方法反演得到的模型剖面；红色虚点线表示新方法反演获得的模型剖面
Fig. 10. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a), and comparison of the S-wave velocity profiles inverted by classical method and new method and the profile of model D (b)
图 11 新方法反演模型的频散函数表面与实测频散曲线的拟合对比
Fig. 11. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve
图 12 基于模型B模拟的含随机噪声的瑞雷波地震记录及其高分辨率频散能量谱
a.含20%随机噪声的地震记录；b.由图a提取的高分辨率能量谱；c.含40%随机噪声的地震记录；d.由图c提取的高分辨率能量谱；图b和d中白色实心点为根据频散能量极大值提取的瑞雷波频散曲线，该频散数据将作为实测频散曲线进行反演
Fig. 12. Rayleigh wave seismic records with random noise and their high-resolution dispersion energy spectra based on model B
图 13 含不同随机噪声数据反演模型剖面与模型B的对比
图中绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型B；红色实线表示无噪声数据反演得到的模型剖面；洋红色点划线表示含20%随机噪声数据反演获得的模型剖面.蓝色点虚线表示含40%随机噪声数据反演获得的模型剖面
Fig. 13. Comparison of the S-wave velocity profiles inverted by different random noise data and the profile of model B
图 14 不同层数设置时反演模型剖面与模型C的对比
图中绿色虚线表示粒子群优化算法反演时模型参数搜索范围；黑色实线表示真实模型C；红色实线表示设置层数为四层时反演得到的模型剖面；洋红色点划线表示设置层数为五层时反演获得的模型剖面.蓝色点虚线表示设置层数为六层时反演获得的模型剖面
Fig. 14. Comparison of the S-wave velocity profiles with different layer setting and the profile of model C
图 15 某高速公路路基瑞雷波勘探实例
a.野外实测的24道瑞雷波地震记录；b.由图 15a中的多道瑞雷波炮集地震记录提取的f-k域频散能量谱.图中白色实线为提取的f-k域实测瑞雷波频散曲线
Fig. 15. Exploration case of Rayleigh wave of a highway roadbed
图 16 传统方法与新方法反演模型正演模拟的频散曲线与实测频散曲线的对比(a)和传统方法与新方法反演的横波速度剖面与钻孔资料的对比(b)
图a中:黑色实点为从图 15中提取得到的实测频散数据变换到f-v域的频散曲线; 红色虚点线为利用新方法反演模型正演模拟获得的多模式频散曲线；蓝色虚线表示利用传统方法反演模型正演模拟获得的基阶波频散曲线.图b中:绿色虚线为粒子群优化算法反演时模型参数搜索范围；带有圆圈的实线表示钻孔资料剖面；红色虚点线表示新方法反演的横波速度剖面；蓝色虚线表示传统方法反演得到的横波速度剖面
Fig. 16. Comparison of the dispersion curves inverted by classical method and new method and the measured dispersion curve (a), and comparison of the S-wave velocity profiles inverted by classical method and new method and the borehole data (b)
图 17 新方法反演模型频散函数表面与实测频散曲线对比拟合情况
Fig. 17. Comparison of the dispersion function surface of model inverted by new method and the measured dispersion curve
表 1 模型A:三层含低速软弱夹层地质模型参数
Table 1. Model A: a three-layer model with a soft layer trapped between two stiff layers
层序号 VS(m/s) VP(m/s) ρ(g/cm3) h(m) 1 220 437 1.8 6 2 160 285 2.0 3 均匀半空间 400 794 2.1 
下载: 导出CSV
表 2 模型B:两层速度递增型地质模型参数及反演搜索范围
Table 2. Model B: a two-layer model characterized by S-wave velocities increasing with depth and search space in the inversion
层序号 模型参数 搜索范围 VS(m/s) VP(m/s) ρ(g/cm3) h(m) VS(m/s) h(m) 1 150 298 1.8 5 100~300 1~10 均匀半空间 450 802 2.1 ∞ 200~3 000 ∞
下载: 导出CSV
表 3 模型C:四层含低速软夹层地质模型参数及反演搜索范围
Table 3. Model C: a four-layer model with a soft layer trapped between two stiff layers and search space in the inversion
层序号 模型参数 搜索范围 VS(m/s) VP(m/s) ρ(g/cm3) h(m) VS(m/s) h(m) 1 200 490 2.0 5 120~300 2~8 2 160 392 2.0 5 100~300 2~8 3 260 637 2.0 5 100~500 2~8 均匀半空间 380 931 2.0 ∞ 200~800 ∞
下载: 导出CSV
表 4 模型D:六层含高速硬夹层地质模型参数及反演搜索范围
Table 4. Model D: a six-layer model with a stiff layer sandwiched between two soft layers and search space in the inversion
层序号 模型参数 搜索范围 VS(m/s) VP(m/s) ρ(g/cm3) h(m) VS(m/s) h(m) 1 200 490 2.0 2.5 150~400 0.1~5 2 260 637 2.0 2.5 10~400 0.1~5 3 120 294 2.0 2.5 50~600 0.1~5 4 240 588 2.0 3 100~600 0.1~5 5 260 637 2.0 3 150~600 0.1~5 6 370 906 2.0 ∞ 200~600 ∞
下载: 导出CSV
-
[1] Cui, J.W., 2004.An Improved Global Optimization Method and Its Application to the Inversion of Surface Wave Dispersion Curves. Chinese Journal of Geophysics, 47(3):521-527(in Chinese with English abstract). doi: 10.1007/s11803-007-0592-y [2] Dal Moro, G., Pipan, M., Gabrielli, P., 2007.Rayleigh Wave Dispersion Curve Inversion via Genetic Algorithms and Marginal Posterior Probability Density Estimation. Journal of Applied Geophysics, 61(1):39-55.doi:10.10 16/j.jappgeo.2006.04.002 [3] Fan, Y.H., Liu, J.Q., Xiao, B.X., 2002.Fast Vector-Transfer Algorithm for Computation of Rayleigh Wave Dispersion Curves. Journal of Hunan University(Natural Sciences Edition), 29(5):25-30(in Chinese with English abstract). doi: 10.1007%2Fs11770-014-0430-8.pdf [4] Feng, H.J., Zhou, A.G., Yu, J.J., et al., 2016.A Comparative Study on Plum-Rain-Triggered Landslide Susceptibility Assessment Models in West Zhejiang Province. Earth Science, 41(3):403-415 (in Chinese with English abstract). doi: 10.1186/s40677-017-0078-9 [5] Foti, S., Comina, C., Boiero, D., et al., 2009.Non-Uniqueness in Surface-Wave Inversion and Consequences on Seismic Site Response Analyses. Soil Dynamics and Earthquake Engineering, 29(6):982-993.doi: 10.1016/j.soildyn.2008.11.004 [6] Gao, L., Xia, J., Pan, Y., 2014.Misidentification Caused by Leaky Surface Wave in High-Frequency Surface Wave Method. Geophysical Journal International, 199(3):1452-1462.doi: 10.1093/gji/ggu337 [7] Haskell, N.A., 1953.The Dispersion of Surface Waves on Multi-Layered Media. Bulletin of the Seismological Society of America, 43(1):86-103.doi: 10.1029/SP030p0086 [8] He, Y.F., Chen, W.T., Chen, X.F., 2006.Normal Mode Computation by the Generalized Reflection-Transmission Coefficient Method in Planar Layered Half Space. Chinese Journal of Geophysics, 49(4):1074-1081(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTOTAL-DQWX200604019.htm [9] Hu, J.F., Duan, Y.K., Hu, Y.L., et al., 1999.Inversion of Shear-Wave Velocity Structure in Shallow Soil from Rayleigh Waces. Chinese Journal of Geophysics, 42(3):393-400(in Chinese with English abstract). http://manu39.magtech.com.cn/Geophy/EN/Y1999/V42/I03/393 [10] Huang, F.M., Yin, K.L., Zhang, G.R., et al., 2015.Landslide Groundwater Level Time Series Prediction Based on Phase Space Reconstruction and Wavelet Analysis-Support Vector Machine Optimized by PSO Algorithm. Earth Science, 40(7):1254-1265(in Chinese with English abstract). http://d.wanfangdata.com.cn/Periodical_dqkx201507013.aspx [11] Huang, Z.X., Li, H.Y., Xu, Y., 2014.Lithospheric S-wave Velocity Structure of West China and Neighboring Areas from Surface Wave Tomography. Chinese Journal of Geophysics, 57(12):3994-4004(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTotal-DQWX201412013.htm [12] Li, Q.C., Shao, G.Z., Liu, J.L., et al., 2006.Past, Present and Future of Rayleigh Surface Wave Exploration. Journal of Earth Sciences and Environment, 28(3):74-77(in Chinese with English abstract). doi: 10.1007/s10712-011-9134-2 [13] Li, X.Y., Zhu, P.M., Zhou, Q., et al., 2014.S-Wave Velocity Structure of Upper Crust in Three Gorges Reservoir Region of the Yangtze River. Earth Science, 39(12):1842-1850(in Chinese with English abstract). http://adsabs.harvard.edu/abs/2014AGUFM.S23C4516L [14] Lin, C.P., Chang, C.C., Chang, T.S., 2004.The Use of MASW Method in the Assessment of Soil Liquefaction Potential. Soil Dynamics and Earthquake Engineering, 24(9-10):689-698.doi: 10.1016/j.soildyn.2004.06.012 [15] Liu, X.F., Fan, Y.H., 2012.On the Characteristics of High-Frequency Rayleigh Waves in Stratified Half-Space. Geophysical Journal International, 190(2):1041-1057.doi: 10.1111/j.1365-246x.2012.05479.x [16] Liu, X.F., Fan, Y.H., Chen, X.F., 2009.Research on the Cross of the Dispersion Curves of Rayleigh Waves and Multi-Modes Coupling Phenomnon. Chinese Journal of Geophysics, 52(9):2302-2309(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTOTAL-DQWX200909016.htm [17] Liu, X.M., Fan, Y.H., Zhai, J.Y., et al., 2009.Numerical Simulation of Rayleigh Wave Zigzag Dispersion Curves of Four Typical Strata. Chinese Journal of Geophysics, 52(12):3042-3050(in Chinese with English abstract). http://www.oalib.com/paper/1569103 [18] Lu, L.Y., Zhang, B.X., Wang, C.H., 2006.Experiment and Inversion Studies on Rayleigh Wave Considering Higher Modes. Chinese Journal of Geophysics, 49(4):1082-1091(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTOTAL-DQWX200604020.htm [19] Pan, J.T., Li, Y.H., Wu, Q.J., et al., 2014.3-D S-wave Velocity Structure of Crust and Upper-Mantle beneath the Northeast China. Chinese Journal of Geophysics, 57(7):2077-2087(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTotal-DQWX201407005.htm [20] Pei, D.H., Louie, J.N., Pullammanappallil, S.K., 2007.Application of Simulated Annealing Inversion on High-Frequency Fundamental-Mode Rayleigh Wave Dispersion Curves. Geophysics, 72(5):R77-R85.doi: 10.1190/1.2752529 [21] Poli, R., Kennedy, J., Blackwell, T., 2007.Particle Swarm Optimization. Swarm Intelligence, 1(1):33-57.doi: 10.1007/s11721-007-0002-0 [22] Shi, Y.L., Jin, W., 1995.Genetic Algorithms Inversion of Lithospheric Structure from Surface Wave Dispersion. Chinese Journal of Geophysics, 38 (2):189-198(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTOTAL-DQWX502.006.htm [23] Shirazi, H., Abdallah, I., Nazarian, S., 2009.Developing Artificial Neural Network Models to Automate Spectral Analysis of Surface Wave Method in Pavements. Journal of Materials in Civil Engineering, 21(12):722-729.doi: 10.1061/(asce)0899-1561(2009)21:12(722) [24] Song, X.H., Tang, L., Lv, X.C., et al., 2012.Application of Particle Swarm Optimization to Interpret Rayleigh Wave Dispersion Curves. Journal of Applied Geophysics, 84:1-13.doi: 10.1016/j.jappgeo.2012.05.011 [25] Tan, K., Zhang, Q.Q., Cao, Q., et al., 2015.Hyperspectral Retrieval Model of Soil Organic Matter Content Based On Particle Swarm Optimization-Support Vector Machines. Earth Science, 40(8):1339-1345(in Chinese with English abstract). http://www.mdpi.com/2072-4292/9/1/28/html [26] Tillmann, A., 2005.An Unsupervised Wavelet Transform Method for Simultaneous Inversion of Multimode Surface Waves. Journal of Environmental & Engineering Geophysics, 10(3):287-294.doi: 10.2113/jeeg10.3.287 [27] Xia, J.H., Miller, R.D., Park, C.B., 1999.Estimation of Near-Surface Shear-Wave Velocity by Inversion of Rayleigh Waves. Geophysics, 64(3):691-700.doi: 10.1190/1.1444578 [28] Xia, J., Chen, C., Li, P.H., et al., 2004.Delineation of a Collapse Feature in a Noisy Environment Using a Multichannel Surface Wave Technique. Géotechnique, 54(1):17-27.doi: 10.1680/geot.54.1.17.36326 [29] Zhang, B.X., Lu, L.Y., Bao, G.S., 2002.A Study on Zigzag Dispersion Curves in Rayleigh Wave Exploration. Chinese Journal of Geophysics, 45(2):263-274(in Chinese with English abstract). http://manu39.magtech.com.cn/Geophy/EN/Y2002/V45/I02/263 [30] Zhang, B.X., Xiao, B.X., Yang, W.J., et al., 2000.Mechanism of Zigzag Dispersion Curves in Rayleigh Exploration and Its Inversion Study. Chinese Journal of Geophysics, 43(4):557-567(in Chinese with English abstract). http://en.cnki.com.cn/Article_en/CJFDTOTAL-DQWX200004016.htm [31] Zhang, S.X., Chan, L.S., 2003.Possible Effects of Misidentified Mode Number on Rayleigh Wave Inversion. Journal of Applied Geophysics, 53(1):17-29.doi: 10.1016/s0926-9851(03)00014-4 [32] Zhou, Z.S., Liu, X.L., Xiong, X.Y., 2007.A Study of Richardson Number and Instability in Moist Saturated Flow. Chinese Journal of Geophysics, 50(2):567-573(in Chinese with English abstract). http://d.wanfangdata.com.cn/Periodical_dqwlxb200702007.aspx [33] 崔建文, 2004.一种改进的全局优化算法及其在面波频散曲线反演中的应用.地球物理学报, 47(3):521-527. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200403025.htm [34] 凡友华, 刘家琦, 肖柏勋, 2002.计算瑞利波频散曲线的快速矢量传递算法.湖南大学学报(自然科学版), 29(5):25-30. http://www.cnki.com.cn/Article/CJFDTOTAL-HNDX200205005.htm [35] 冯杭建, 周爱国, 俞剑君, 等, 2016.浙西梅雨滑坡易发性评价模型对比.地球科学, 41(3):403-415. http://www.earth-science.net/WebPage/Article.aspx?id=3259 [36] 何耀锋, 陈蔚天, 陈晓非, 2006.利用广义反射-透射系数方法求解含低速层水平层状介质模型中面波频散曲线问题.地球物理学报, 49(4):1074-1081. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200604019.htm [37] 胡家富, 段永康, 胡毅力, 等, 1999.利用Rayleigh波反演浅土层的剪切波速度结构.地球物理学报, 42(3):393-400. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX199903011.htm [38] 黄发明, 殷坤龙, 张桂荣, 等, 2015.基于相空间重构和小波分析-粒子群向量机的滑坡地下水位预测.地球科学, 40(7):1254-1265. http://www.earth-science.net/WebPage/Article.aspx?id=3113 [39] 黄忠贤, 李红谊, 胥颐, 2014.中国西部及邻区岩石圈S波速度结构面波层析成像.地球物理学报, 57(12):3994-4004. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX201412013.htm [40] 李庆春, 邵广周, 刘金兰, 等, 2006.瑞雷面波勘探的过去、现在和未来.地球科学与环境学报, 28(3):74-77. http://www.cnki.com.cn/Article/CJFDTOTAL-XAGX200603016.htm [41] 李小勇, 朱培民, 周强, 等, 2014.三峡库区上地壳横波速度结构.地球科学, 39(12):1842-1850. http://www.earth-science.net/WebPage/Article.aspx?id=3000 [42] 刘雪峰, 凡友华, 陈晓非, 2009.Rayleigh波频散曲线"交叉"及多模式耦合作用研究.地球物理学报, 52(9):2302-2309. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200909016.htm [43] 刘雪明, 凡友华, 翟佳羽, 等, 2009.四种典型地层的瑞雷波"之"字型频散曲线数值模拟研究.地球物理学报, 52(12):3042-3050. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200912014.htm [44] 鲁来玉, 张碧星, 汪承灏, 2006.基于瑞利波高阶模式反演的实验研究.地球物理学报, 49(4):1082-1091. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200604020.htm [45] 潘佳铁, 李永华, 吴庆举, 等, 2014.中国东北地区地壳上地幔三维S波速度结构.地球物理学报, 57(7):2077-2087. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX201407005.htm [46] 石耀霖, 金文, 1995.面波频散反演地球内部构造的遗传算法.地球物理学报, 38(2):189-198. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX502.006.htm [47] 谭琨, 张倩倩, 曹茜, 等, 2015.基于粒子群优化支持向量机的矿区土壤有机质含量高光谱反演.地球科学, 40(8):1339-1345. http://www.cnki.com.cn/Article/CJFDTOTAL-DQKX201508009.htm [48] 张碧星, 鲁来玉, 鲍光淑, 2002.瑞利波勘探中"之"字形频散曲线研究.地球物理学报, 45(2):263-274. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200202012.htm [49] 张碧星, 肖柏勋, 杨文杰, 等, 2000.瑞利波勘探中"之"形频散曲线的形成机理及反演研究.地球物理学报, 43(4):557-567. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200004016.htm [50] 周竹生, 刘喜亮, 熊孝雨, 2007.弹性介质中瑞雷面波有限差分法正演模拟.地球物理学报, 50(2):567-573. http://www.cnki.com.cn/Article/CJFDTOTAL-DQWX200702029.htm 期刊类型引用(8)
1. 宫丰，陈晓非，凡友华，刘雪峰，唐海兵. 基于局部优化粒子群算法的背景噪声反演. 地球物理学报. 2023(06): 2550-2563 . 
百度学术2. 李巧灵，雷晓东，朱琳，李晨，陈蓓蓓. 基于微动探测地下空间异常响应特征研究. 工程勘察. 2022(01): 70-78 . 百度学术3. 崔岩，王彦飞. 瑞雷波多阶模式频散曲线稀疏正则化反演方法研究. 地球物理学报. 2022(03): 1086-1095 . 百度学术4. 宫丰，陈晓非，凡友华，刘雪峰，唐海兵. 久期函数值作为目标函数的瑞雷波频散曲线反演. 工程地球物理学报. 2022(05): 742-752 . 百度学术5. 王序宇，张立，孙康，陈雨翔，冯栋，陈善林，蒋吉佑，卢柱峰，刘争平. 水平层状介质中Rayleigh波多模式速度频散和椭圆极化频散特征物理实验模拟研究. 地球物理学进展. 2022(05): 2044-2052 . 百度学术6. 董智开，段文胜，肖承文，胡天跃，张献兵. 基于快速标量传递算法的瑞雷波频散曲线反演研究. 北京大学学报(自然科学版). 2020(04): 614-628 . 百度学术7. 何达喜，朱培民，申随水，蔡伟，汪文刚. 煤层物性参数对瑞雷型槽波频散曲线“交叉”现象的影响. 煤矿安全. 2020(12): 69-74 . 百度学术8. 王书纯，黄用勤，葛健，陈珺. 高精度分布式面波地震仪. 中国测试. 2018(06): 66-70 . 百度学术其他类型引用(5)
-
下载:
百度学术
点击查看大图
计量
- 文章访问数: 7859
- HTML全文浏览量: 1702
- PDF下载量: 36
- 被引次数: 13
