地球物理学报  2011, Vol. 54 Issue (11): 2951-2959   PDF    
基于改进粒子群算法的地震标量波方程反演
朱童1,2, 李小凡1, 李一琼1,3, 张美根1     
1. 中国科学院地球深部重点实验室,中国科学院地质与地球物理研究所,北京 100029;
2. 中国石油化工股份有限公司石油物探技术研究院,南京 210014;
3. 中国地震局地球物理研究所,北京 100081
摘要: 针对标准粒子群优化(PSO)算法存在易出现早熟而陷入局部最优以及进化后期收敛速度慢等缺陷,通过考虑粒子所处位置间相互作用,提出了一种改进的并行粒子群优化算法.由于引入粒子位置间的相互影响,减少了粒子搜索过程盲目性,因此能有效提高算法的收敛速度.数值试验表明,这种改进的粒子群算法适用于二维标量波方程的速度反演,且算法具有对初始模型依赖性低、收敛速度快、反演结果稳定、抗噪能力强等特点,为进一步将该反演算法用于弹性波波动方程以及弹性参数反演提供了理论依据.
关键词: 改进粒子群算法      标量波方程      速度反演      并行     
Seismic scalar wave equation inversion based on an improved particle swarm optimization algorithm
ZHU Tong1,2, LI Xiao-Fan1, LI Yi-Qiong1,3, ZHANG Mei-Gen1     
1. Key Laboratory of the Earth's Deep Interior,CAS, Institute of Geology and Geophysics, Chinese Academy of Sciences, Beijing 100029,China;
2. SINOPEC Geophysical Research Institute, Nanjing 210014, China;
3. Institute of Geophysics, China Earthquake Administration, Beijing 100081,China
Abstract: In the standard particle swarm optimization (PSO), the premature convergence of particles and slow convergence in the late process decrease the searching ability of the algorithm. In this paper, by taking the positions of the particles into consideration, we propose an improved parallel particle swarm optimization (IPPSO) algorithm, which can increase the efficiency, to reduce the blindness in the search process. The performance of this improved particle swarm optimization method in solving 2-D scalar wave equation inversion problems is investigated. Our numerical experiments indicate that this method is suitable for scalar velocity inversion problems since it has merits such as low dependence on initial model, high constringency speed, stable result and strong antinoise ability. The results hold promise for further elastic wave equation and elastic parameters inversion studies.
Key words: Improved particle swarm optimization      Scalar wave equation      Velocity inversion      Parallel     
1 引言

波动方程反演问题的研究具有重要的理论意义和实用价值,在地球物理学等领域有着重要应用.然而,由于波动方程反演问题的目标函数是一个高度非线性的非凸复杂多峰函数,受到非线性与不适定性的双重困扰,加之实际的物理现象的复杂性,使得波动方程反演成为地球物理反演领域的一个热点和难点问题.反演过程具有不适定性、局部收敛性、计算量大、高阶导数难于求解、初始猜测很难合理选择的特点,所以反演方法的优劣就成了反演成功与否的关键.因此,开展实用、可靠的数值反演方法的研究具有十分重要的意义.

自20世纪70年代以来,许多学者从不同角度探讨波动方程的反演,提出了许多有特色的方法,这些方法都有各自的优缺点:脉冲谱方法(Pulse-Spectrum Technique, PST)是由Tsien和Chen[1]于1974年在求解流体力学波速反演问题时提出的一种数值反演方法,基本思想是在一个具紧支集的时域中对测量数据采样,而在复频域中用迭代法求出未知参数的数值解.由它发展起来的广义脉冲谱方法(GPST)[2]在一定具体形式下是收敛的,具有对原始信息需求量较小等优点;不足之处在于GPST实质上是牛顿型迭代法,因此是局部收敛的;正则化参数选取主要凭经验,往往很难达到最优迭代;另外,GPST 需要反复求解正问题,计算量也不小[3].Born反演方法是由Cohen和Bleistein[4]在1977年提出的一种求解逆散射问题的扰动方法,随后Tarantola[5]、Ikelle 等[6]和Zhang[7]基于常背景的Born近似,给出了波动方程反演的一些方法,以及在此基础上发展起来的Born-WKBJ[8]方法.该类方法将一个非线性反演问题转化为线性问题,使反演过程大大简化,同时计算量也大大减少;其缺点是要求给出一个非常接近真实地层的初始模型,此外要求扰动参数必须很小,当扰动参数较大时,反演结果误差较大.而只有当初始的速度模型在目标函数的全局极小点的邻域内时,近似线性化的关系才能成立.针对波动方程的奇异性沿特征线传播的特点,Xie等[9]提出了时卷特征正则化方法(TccR),应用到散射势反问题与P 波速、声波系数反演问题当中.此方法的主要思想是构造恰当的卷积算子,将原来的微分方程反问题转化为非线性积分几何方程,再设法推得其FreChet微分算子来构造Newton迭代格式.数值模拟结果表明此法比较稳定可靠,而且计算量大大减少,但它存在迭代方法所固有的通病即局部收敛性.一些学者采用基于最小平方准则的完全非线性反演方法[1011],先给出初始模型,在泛函空间中求正问题的解与观测数据之间的拟合差作为目标函数,将波动方程反问题转化为非线性优化问题,然后使用非线性最优化算法包括共轭梯度法、最速下降法[12]、拟牛顿型迭代法[13]、伴随状态法[14],通过不断地修改模型参数,直到收敛到使得目标函数达到最小的参数值,从原理上更适合于含有多次反射的观测数据.反问题的非线性反演算法的最大优点是通用性强,不受维数的限制,有一定的抗噪声能力.其缺点是由于所采取的极小化准则的限制,获得的解不能保证是整体最优;一般要求初始模型比较接近真实模型,另外该方法的计算量也是非常大的.Rotmhna和Smiht等[15]分别将统计学中的模拟退火方法(Simulated Annealing)[16]和多参数优化中的遗传算法(GeneticAlgorithms)应用到非线性反演中用来克服用一般的优化方法所进行的非线性反演的局部极小值问题,但是也难以突破计算量大和易入局部极值的瓶颈.冯国锋等[17]发展了一种波动方程问题的多尺度信赖域的反演方法,用信赖域的思想来解决反演问题,但对于方法的收敛性和收敛速度在理论上未作分析.

粒子群优化算法(Particle Swarm Optimization, PSO)[18]是一种基于群体的随机全局优化算法,最早由美国社会心理学家R.Eberhart和J.Kennedy于1995年提出.与其他基于群体的进化算法不同的是,它模拟的是社会行为:将每个可能产生的解表述为群体中的一个微粒,并为之规定简单的行为规则,从而使整个粒子群表现出复杂特性,用来解决复杂的优化问题.由于该优化算法易于理解实现,可以解决非线性问题,在短期内得到很大发展,已在函数寻优、神经网络训练、数字电路优化、TSP 问题以及工程应用等多个优化及反演领域得到了成功的应用,但目前在地球物理反演领域[1920]的应用尚不多见,本文尝试将粒子群算法用于解决地震标量波方程反演这一高度非线性多极值优化问题.

2 粒子群算法 2.1 基本PSO算法

PSO 是一种基于群体的优化算法,该算法源于对鸟群捕食行为的研究,通过对简单社会系统的模拟,在多维解空间中构造所谓的“粒子群",粒子群中每个粒子通过跟踪自己和群体所发现的最优值,修正自己的前进方向和速度,从而实现寻优.群体中每个粒子表示问题的一个可行解,并具有与目标函数相关的适应度值.粒子在搜索空间中以一定的速度飞行,并根据自身的飞行经验以及当前最优粒子的状态对速度进行动态调整,个体之间通过协作与竞争,实现对问题最优解的搜索.

设在n维搜索空间中,由m个粒子组成种群X= {x1,…,xi,…,xm},其中第i个粒子位置为xi= (xi1xi2,…,xin)T,速度为vi= (vi1vi2,…,vin)T.粒子i经历的最好位置记作pi= (pi1pi2,…,pin)T,整个粒子群经历的最好位置记作pg =(pg1pg2,…,pgn)T,粒子xi将按式(1)、式(2)改变速度和位置.

(1)

(2)

粒子i通过式(1)、式(2)决定下一步的运动位置.式中,s=1,2,…,n,代表搜索空间的维数,ω 为惯性权值,代表前一时刻的速度对当前时刻速度的影响.ω 较小时,粒子容易陷入局部最优,当ω 较大时,粒子可以跳出局部极值,但收敛速度较慢.因此调整ω 的大小可以在粒子收敛性和计算速度之间取得平衡,通常在迭代初期ω 取值较大,收敛速度较快,ω 随着迭代过程线性变小,在迭代后期使得解跳出局部极值.c1c2 为系数,c1 为粒子自身的加速度权重,代表粒子向自己经验的学习,c2 为全局加速度权重,代表粒子向群体经验的学习;r1r2 为介于[0, 1]间的相互独立的随机数,一般服从正态分布,从而动态调节粒子向自身经验学习和向群体经验学习的权重.为了减少进化过程中粒子离开搜索空间的可能,通常对速度vi限定在一个范围[vmin, vmax]之内,而粒子的位置通常由先验信息加以约束.

2.2 粒子群算法的改进

粒子群算法是全局搜索算法,实现简单,然而基本粒子群算法也存在着收敛速度慢、容易陷入局部极值等缺点,许多学者都对基本粒子群算法进行过改进,包括算法本身的改进[2122]如更新新的速度位置公式以及与其他进化方法的结合等[2324].

本文中,考虑粒子前一次位置对更新之后得到的新位置的影响,引入了影响因子α,对基本粒子群算法公式进行了更新.即在公式(1)、(2)后通过

(3)

来进行粒子新位置的选择.

式(3)中,xis(t+1)为粒子的当前位置,xis(t)为粒子的前一次位置,xis′(t+1)为粒子的新位置.α介于[0, 1]之间,代表粒子的前一次位置和当前位置在粒子新位置选择中所占的权重.

α 较小时,粒子的前一次位置对新位置的影响较小,算法的收敛速度和基本粒子群算法的收敛速度相当;当α 逐渐增大时,粒子的前一次位置对新位置的选择所占的权重加大,使得粒子在选择新位置的过程中有了参考点,减少了粒子位置选择的盲目性,因此收敛速度会逐渐提高;而当α 过大时,粒子由于受到前一次位置的影响过大,则往往容易陷入局部极值.因此,算法中,α 的取值是影响算法收敛速度的关键因素.数值试验验证了上述分析过程.数值试验的结果表明[25]:当α≤10%时,对算法收敛速度提高的不明显;当α≥50%时,算法容易陷入局部极值;而当10%<α<50%时,算法的收敛速度基本随着α 的增大而提高.

事实上,如果我们将公式(2)带入公式(3),可以得到

(4)

这表明,对位置的影响因子α,相当于对速度的约束因子(1-α).这相当于在搜索空间里面对搜索的方向未加限制,而对搜索的步长加以约束.这样能够有效地避免在搜索过程中的由于步长过大而错过了全局极值,减少了搜索过程中的盲目性,因此能提高算法的搜索效率.而文献[26]也从数学角度证明了通过对搜索过程中的速度约束,能保证算法的收敛性和提高搜索效率.本文与文献[26]区别在于从不同角度提出的对速度约束的策略,而且约束因子的取值范围不同.相同测试函数结果表明本文提出的改进算法的约束因子相比文献中的算法的约束因子的计算效率更高,因此本文将这种改进的粒子群算法用于波动方程反演问题.

3 波动方程反演问题 3.1 数学模型

描述地震波传播的二维声波方程[27]

(5)

其中u(xzt)为位移函数,v(xz)为介质在(xz)点的速度.若考虑平面波入射,方程变为

(6)

其中,G(t)为震源函数,δ(z)为δ 函数.在区域Ω =[0,X]× [0,Z]上考虑方程,边界条件可描述为

(7)

初始条件为

(8)

(6) ~(8)式组成了声波方程正演模型,即由已知速度v(xz)可解出u(x,0,t)=f(xt),其中f(xt)相当于地表的地震记录,反演的目的是已知f(xt)来求介质速度v(xz).

对式(5)、(6)进行差分离散化,得到速度和地震记录的离散量vijfik.按照适当的顺序排列,组成向量的形式VF.令

(9)

地震记录组成了

(10)

于是反演vij问题变成求如下泛函的问题:

(11)

3.2 波动方程的粒子群反演

对于波动方程反演问题,通过定义目标函数

(12)

可将式(11)的泛函反演问题转化为求解优化问题.式中V= (v1v2,…,vM)为介质参数向量(如标量波速度),M为总的参数个数(离散后的网格点总数);FiFi′(i=1,2,…,N)分别为观测数据和理论模型模拟的地震记录,N表示按一定规则排列后地震记录向量的数量.取目标函数式(12)用来作为模型拟合度函数时,问题可表述为:求σ*X,使

(13)

其中,

3.3 波动方程的粒子群反演步骤

本文中,将待反演的介质速度对应粒子群算法中粒子的“位置",将介质速度的扰动变化对应为算法中的粒子的“速度".给定初始模型后和先验信息限定解的变化范围后,通过粒子群算法的规则,从而不断修正初始模型,最终得到问题的解.

类似基本粒子群算法的实现步骤,可以得到粒子群算法解决波动方程问题的步骤为

① 设置算法的参数和最大迭代次数N:本文中参数为ωmax=0.9,ωmin=0.4,c1=c2=2,α=0.4,最大迭代次数Nmax为2000,ω 变化截止迭代次数为1500(当迭代次数N从1增大至1500 时,ω 由0.9线性减小为0.4,N>1500时,ω=ωmin=0.4);

② 在满足控制变量约束条件下赋予粒子群中每个粒子初始速度和初始位置;

③ 对群体中的每个粒子,利用式(5)、(6)正演出地震记录,然后通过式(12)计算它们的适应度值,找到每个粒子经历的最好位置Pi 和群体经历的最好的位置为Pg;

④ 利用式(1)、(2)和(3)更新粒子速度和位置;

⑤根据先验信息对粒子位置和速度进行约束处理:若大于最大值,赋为最大值;若小于最小值,赋为最小值;

⑥ 判断是否满足迭代收敛条件,若不满足结束条件,则n=n+1,并执行步骤③;若满足结束条件,则执行步骤⑦;

⑦ 输出最优位置Pg.

最优位置Pg 即对应为反演出的介质速度.

3.4 并行粒子群算法的实现

在实际计算中,由于波动方程反演问题的正演计算量较大,单机难以完成,而粒子群算法又具有天然的并行性,加上并行编程环境MPI可以提供组通信使各个进程共享规约结果[28],实现起来较为简单,因此本文采用了MPI并行编程环境实现了改进的粒子群算法的并行化,把全局的粒子个数分配到各个进程进行处理,然后通过组规约找出子迭代过程中的最优解以及所在的位置记为当前迭代次数下的全局最优解,通过组广播,以全局最优解更新各子群最优解.改进后的并行粒子群算法流程图如图 1所示.其中,子群内的改进粒子群算法迭代的步骤如3.3节所述.由于粒子种群的多样性,将粒子分配到各个子进程进行计算并不会增加额外的计算量.而算法调用了MPI组通信的规约操作,各个进程能够共享规约结果,通过进程阻塞来实现同步.因此并行的通讯仅集中在各子群搜索到子最优解后的规约广播和阻塞操作,通讯开销极小,而且能保证粒子种群的多样性.当某一进程搜索到更接近的全局最优解时,会带领其他进程一起向全局最优解靠近搜索,从而大大提高算法的搜索效率.

图 1 改进的并行粒子群算法流程图 Fig. 1 Flow diagram of improved parallel PSO algorithm
4 数值试验

本节中,通过几个数值试验来验证上述算法用于标量波方程反演的正确性和有效性.

4.1 水平层状模型

首先考虑水平层状介质,设计如图 2所示的模型1.hx=hz=10m, x方向、z方向均取64个网格,时间采样间隔为0.001s, 时间采样为300ms.震源为主频50Hz的雷克子波,分别位于(16,32),(32,32),(48,32)处.检波器位于地表处水平接收.三层速度分别为2700m·s-1、3000m·s-1和2500m·s-1,初始速度模型为2500 m·s-1的均匀模型.速度取值约束范围为[2500, 3000],速度变化扰动范围为[-200, 200].粒子群算法的参数与3.3节中所述一致,每个子群由8个粒子组成.正演采用带三阶时间保辛结构的伪谱法[29],采用浪潮TS1000高性能并行机组,16台计算节点参与反演计算.

图 2 模型1三层水平均匀模型 Fig. 2 Model 1 3-layer homogeneous model

反演结果如图 3 所示.从反演结果可以看出,该粒子群算法适用于地震标量波方程的速度反演.速度分界层清晰可见,反演得到的各层速度值也十分逼近真值.由于采用并行算法,反演效率较高,约耗时2h.

图 3 模型1的反演结果 Fig. 3 Inversion resutt of model 1

图 4为反演得到的模型的地震记录,图 5为反演模型与真实模型的地震记录残差图.图 6 为检波器在地表(1,32)处接收到的地震记录,图 5、6 同样说明了反演得到的模型和真实模型很接近.

图 4 反演得到的模型的地震记录 Fig. 4 Seismogram of inversion model
图 5 反演模型与真实模型的地震记录残差图 Fig. 5 Residual between mversion model and true model
图 6 地表处的地震记录 实线表示真实模型,“ + ”线表示反演模型,虚线表示残差 Fig. 6 Seismogram in surface line: true model,dash line with mark “ + ”: inversion model,dash line: residual

为了研究不同初始模型对反演结果的影响,设计初始模型分别为3000m·s-1的均匀模型和一个速度值介于2500m·s-1和3000 m·s-1之间的随机初始模型,两个模型的反演结果如图 7、8 所示.从反演结果来看,两个不同的初始模型得到的反演结果与前面得到的结果(图 3)相差无几,速度分界面和速度值均能够较好地反演出.图 3图 7图 8可以说明,初始模型的选择对反演结果的影响不大,说明该算法对初始模型的依赖性极小.

图 7 初始模型为3000 m • s-1的反演结果 Fig. 7 Inversion result of the initial model (3000 m • s-1)
图 8 初始模型为随机模型的反演结果 Fig. 8 Inversion result of random initial model

为了进一步研究该算法的抗噪能力,我们对模型1生成的地震记录添加了10%的随机噪声,然后进行反演.图 9为反演后的结果.可以看出,速度分界层清晰可见,除第一层反演得到的速度有一定的波动外,第二层介质和第三层介质得到的结果与真实值十分接近.这是由于反演过程中受到了添加噪声的影响,而下面两层介质由于真实速度分别是反演过程的速度取值范围的极大值和极小值,因此噪声对其影响并不明显.分析第一层速度数据可以得出,噪声造成的最大速度扰动也仅为背景速度的1%左右,这说明该算法具有良好的抗噪能力.

图 9 模型1加上10%噪声的反演结果 Fig. 9 Inversion result of model 1 with 10% noise
4.2 横向非均匀模型

为了验证算法对横向非均匀模型的反演能力,设计如图 10所示的模型2,上层为3000m·s-1,下层速度为2500 m·s-1,速度分界层为倾斜粗糙界面,初始模型为2500m·s-1和3000m·s-1的均匀模型.其他参数与模型1一致.由于此模型为横向非均匀模型,3.4 节中介绍的全局广播的并行方法并不完全适用.因此在实际反演过程中,对整个反演区域进行分区处理,除了通过全局规约、广播找到全局最优解外,每个子进程还与周围邻近的进程进行通讯,通过比较与周围进程最优解的差值和与全局最优解的差值大小来判断是取全局最优还是邻近的局部最优.用通过这样处理后的并行粒子群算法进行了反演,图 11a图 11b是分别是初始模型为2500m·s-1和3000m·s-1对应的反演结果.由于数据离散时造成的粗糙界面会产生绕射、散射等复杂波形,这在一定程度上会加大反演的难度,加上真实模型生成的时候是连续的,而反演的时候由于对区域进行了分块处理,与真实模型的界面不能完全匹配,所以反演出的分界面与真实模型并未完全重合,但分层界面位置正确且清晰可见,分界层两边速度对比明显,反演的速度也与真实模型较为接近,说明该反演方法也适用于非均匀模型的反演.两个不同的初始模型均能反演出较为理想的结果,说明该方法对初始模型依赖性很小.

图 10 模型2两层横向非均匀模型 Fig. 10 Model 2 2-layer lateral inhomogeneous model
图 11 模型2的反演结果 (a)初始模型为2500 m • s-1的均匀模型;(b)初始模型为3000 m • s-1的均匀模型. Fig. 11 Inversion results of model 2 (a) Initial homogenous model, the velocity is 2500 m • s-1 ; (b) Initial homogenous model, the velocity is 3000 m • s-1

设计如图 12所示的三层非均匀模型3,模型含一U 型凹陷粗糙分界面.三层速度分别为3000m·s-1,2700m·s-1和2500m·s-1.初始模型为2500m·s-1的均匀模型.其他参数如前所述.

图 12 模型3三层非均勻模型 Fig. 12 Model 3 3-layer lateral inhomogeneous model

图 13为该模型的反演结果.从图上可以看出,对于界面分界层和大部分区域的速度值反演的结果都很理想,反演结果不理想的区域集中在U 型凹陷分界面附近及下方.这是由于网格离散时造成的倾斜界面并不是光滑的,地震波传播时遇到粗糙分界面以及界面的拐角点时会产生散射波等复杂波形,与透射波形交织在一起,使得地表附近的检波器接收到的波形数据也变得复杂,增加了反演的难度.加之实际反演过程中分区域处理,使得反演非均匀介质时很难与真实模型完全吻合.尽管如此,仍然能从图中看出U 型分解面的轮廓以及分界层两边的速度值对比.说明该算法同样适用于较为复杂的模型反演.图 13图 11b分别对应在不同初始模型条件下反演不同非均匀模型,它们间的对比也进一步验证了在反演非均匀模型介质时,该方法对初始模型依赖性小.

图 13 模型3反演结果 Fig. 13 Inversion resutt of model 3
5 结论和讨论

本文对基本粒子群算法进行了改进,通过考虑粒子更新时位置的相互影响,提出了一种改进的粒子群算法,并将其并行化以用于二维地震标量波波动方程速度反演.由于采用了并行群体搜索策略,它反演收敛的速度快;不依赖于初始模型,无论初始模型是低速层还是高速层,甚至是随机速度模型,均可以反演出理想的结果;且算法具有较强的抗噪能力.通过对水平层状模型和横向非均匀模型进行反演,结果表明,本文构造的这种改进的粒子群算法用于标量波方程速度反演是可行的,且算法具有对初始模型依赖性小、收敛速度快、计算效率高等优点,为进一步将该反演算法用于弹性波波动方程以及弹性参数反演提供了理论依据.

取仅为通过前人文献以及数值试验试算取得,并未在理论上加以过多分析算法参数的选取对算法收敛性的影响.此外,对于非均匀复杂模型的反演,文中所设计的反演方法有待于进一步优化.后续工作中,我们将进一步分析粒子群算法参数的选取、较大规模的复杂模型反演、非均匀模型反演方法优化以及多参数、弹性波波动方程的反演等问题.

参考文献
[1] Tsien D S, Chen Y M. A numerical method for nonlinear inverse problems in fluid dynamics. Proc. Int. Conf. Comput. Meth. Nonlinear Mechs. Austin: Univ of Taxas Press, 1974 : 935 -943.
[2] Chen Y M. Generalized pluse-spectrum technique. Geophysics , 1985, 50(11): 1664-1675. DOI:10.1190/1.1441857
[3] Han B, Feng G F, Liu J Q. A widely convergent generalized pulse-spectrum technique for the inversion of two-dimensional acoustic wave equation. Applied Mathematics and Computation , 2006, 172(1): 406-420. DOI:10.1016/j.amc.2005.02.011
[4] Cohen J K, Bleistein N. An inverse method for determining small variations in propagation speed. SIAM J. Appl. Math. , 1977, 32(4): 784-799. DOI:10.1137/0132066
[5] Tarantola A. Linearized inversion of seismic reflection data. Geophysical Prospecting , 1984, 32(6): 998-1015. DOI:10.1111/gpr.1984.32.issue-6
[6] Ikelle L T, Diet J P, Tarantola A. Linearized inversion of multi-offset seismic reflection data in the w-k domain. Geophysics , 1986, 51(6): 1266-1276. DOI:10.1190/1.1442179
[7] Zhang X K. Generalized inversion and its application in inverse scatting problems. J. Comput. Math. , 1989, 7(4): 374-382.
[8] Clayton R W, Stolt R H. A Born-WKBJ inversion method for acoustic reflection data. Geophysics , 1981, 46(11): 1559-1567. DOI:10.1190/1.1441162
[9] Xie G Q. A new iterative method for solving the coefficient inverse problem of the wave equation. Comm. Pure and Appl. Math. , 1986, 39(3): 307-322. DOI:10.1002/(ISSN)1097-0312
[10] Mara P. Nonlinear two-dimensional elastic inversion of seismic reflection data. Geophysics , 1987, 52(9): 1211-1228. DOI:10.1190/1.1442384
[11] Tarantola A. Theoretical background for the inversion of seismic waveforms, including elasticity and attenuation. Pure and Applied Geophysics , 1988, 128(1-2): 365-399. DOI:10.1007/BF01772605
[12] Bleistein N, Cohen J K, Hagin F G. Two and one-half dimensional born inversion with an arbitrary reference. Geophysics , 1987, 52(1): 26-36. DOI:10.1190/1.1442238
[13] Xie G Q, Li J H, Chen Y M. Gauss-Newton-Regulating method for solving coefficient inverse problem of partial differential equation and its convergence. J. Comput. Math. , 1987, 5(1): 38-49.
[14] 孙豪志, 范祯祥. 弹性波介质中纵、横波速度反演. 石油地球物理勘探 , 1994, 29(6): 678–684. Sun H Z, Fan Z X. P and S-wave velocity inversions in elastic medium. Oil Geophysical Prospecting (in Chinese) (in Chinese) , 1994, 29(6): 678-684.
[15] Sen M K, Stoffa P L. Nonlinear one-dimensional seismic waveform inversion using simulated annealing. Geophysics , 1991, 56(10): 1624-1638. DOI:10.1190/1.1442973
[16] Stoffa P L, Sen M K. Nonlinear multi-parameter optimization using genetic algorithms: Inversion of plane-wave seismograms. Geophysics , 1991, 56(11): 1794-1810. DOI:10.1190/1.1442992
[17] 冯国锋. 波动方程反问题的多尺度—信赖域反演方法. 哈尔滨: 哈尔滨工业大学 , 2006. Feng G F. Multi-Scale Trust Region Inversion Methods for Inverse Problem of Wave Equations (in Chinese). Harbin: Harbin Institute of Technology (in Chinese) , 2006.
[18] Eberhart R, Kennedy J. A new optimizer using particle swarm theory. In: IEEE. Proceedings of the Sixth International Symposium on Micro Machine and Human Science. Piscataway: IEEE Service Center, 1995. 39~43
[19] 王福昌, 王永革, 胡顺田. 粒子群算法在主震断层面参数估 计中的应用. 地震研究 , 2008, 31(2): 149–153. Wang F C, Wang Y G, Hu S T. Application of particle swarm optimization to the estimation of mainshock fault plane parameters. Journal of Seismological Research (in Chinese) (in Chinese) , 2008, 31(2): 149-153.
[20] 李刚毅, 蔡涵鹏. 基于粒子群优化算法的波阻抗反演研究. 勘探地球物理进展 , 2008, 31(3): 187–191. Li G Y, Cai H P. Wave impedance inversion based on particle swarm optimization. Progress in Exploration Geophysics (in Chinese) (in Chinese) , 2008, 31(3): 187-191.
[21] Shi Y, Eberhart R. A modified particle swarm optimizer. In: Proceeding of the IEEE Conference on Evolutionary Computation. Piscataway. NJ: IEEE Press , 1998: 69-73.
[22] Shi Y, Eberhart R. Fuzzy adaptive particle swarm optimization. In: Proceeding of the IEEE Conference on Evolutionary Computation. Seoul, Korea, 2001. 101~106
[23] Natsuki Higashi, Hitoshi Iba. Particle swarm optimization with Gaussian mutation. In: Proceedings of the Swarm Intelligence Symposium, SIS2003 and IEEE , 2003.
[24] Lvbjerg M, Rasmussen T K, Krink T. Hybrid particle swarm optimiser with breeding and subpopulations. In: Proceedings of the Genetic and Evolutionary Computation Conference , 2001.
[25] 朱童, 李小凡, 鲁明文. 位置加权的改进粒子群算法. 计算机工程与应用 , 2011, 47(5): 4–6. Zhu T, Li X F, Lu M W. Improved particle swarm optimization algorithm with position weighted. Computer Engineering and Applications (in Chinese) (in Chinese) , 2011, 47(5): 4-6.
[26] 张建科. 几类改进的粒子群算法. 西安: 西安电子科技大学 , 2006. Zhang J K. Several Classes of Improved Particle Swarm Optimization (in Chinese). Xi'an: Xidian University (in Chinese) , 2006.
[27] 王守东, 刘家琦, 黄文虎. 二维声波方程速度反演的一种方法. 地球物理学报 , 1995, 38(6): 833–839. Wang S D, Liu J Q, Huang W H. A method of velocity inversion of two dimensional acoustic wave equation. Chinese J. Geophys. (Acta Geophysica Sinica) (in Chinese) (in Chinese) , 1995, 38(6): 833-839.
[28] 赵勇, 岳继光, 李炳宇, 等. 一种新的求解复杂函数优化问题的并行粒子群算法. 计算机工程与应用 , 2005, 41(16): 58–60. Zhao Y, Yue J G, Li B Y, et al. A parallel particle swarm optimization algorithm based on multigroup for solving complex functions optimization. Computer Engineering and Applications (in Chinese) (in Chinese) , 2005, 41(16): 58-60.
[29] 李一琼, 李小凡, 朱童. 基于辛格式奇异核褶积微分算子的地震标量波场模拟. 地球物理学报 , 2011, 54(7): 1827–1834. Li Y Q, Li X F, Zhu T. The seismic scalar wave field modeling by symplectic scheme and singular kernel convolutional differentiator. Chinese J. Geophys. (in Chinese) , 2011, 54(7): 1827-1834.