基于梯度投影法的全变差正则化全波形反演

【www.zhangdahai.com--其他范文】

姚含,徐海

(贵州省地质环境监测院,贵州 贵阳 550001)

我国油气资源勘探从构造勘探阶段逐步走向岩性勘探阶段,勘探难度逐渐加大。地震勘探具有勘探深度大、精度高、数据量大等特点,是目前油气勘探最有效的方法,应用广泛。为进一步提高地震勘探精度,地震全波形反演方法迅速发展起来,成为国内外学者研究的热点。

全波形反演(full waveform inversion,FWI)是一种高精度、高分辨率的速度反演方法[1-2]。FWI利用地震资料中的运动学和动力学信息,通过最小化模拟数据与实际观测数据之间的误差来重构地下介质模型[3-4],能得到地下介质中更多细节的参数变化特征[5]。由于实际观测中,子波未知,数据频带有限,且通常存在一定程度的噪声干扰,波场与介质参数之间呈现较强的非线性关系,FWI体现为高度欠定问题,对初始模型的精度要求较高,导致FWI在实际应用中达不到预期的效果。

为了减少FWI的不适定性,加入先验模型约束的正则化项,能使反演具有更好的稳定性[6-7]。全变差(total variation,TV)约束方法是一种非光滑约束,由Rudin等[8]在1992年提出,主要用于图像的去噪处理,能有效处理不连续性解的求解问题,重构图像中的一些间断点,保留边缘信息。近几年来,全变差约束已经应用在电法成像[9]、图像去噪[10-11]、偏移速度分析[12-13]、混合数据偏移[14]等领域。地下地层通常存在较多复杂的构造变化带及岩性突变带,导致地层物理参数在这些地方产生不连续变化,应用基于TV约束的FWI可以有效保留这些不连续特征[15]。

本文采用梯度投影法[16],构建基于全变差约束的全波形反演目标函数,实现多约束集下的声波方程全波形反演,逐步搜索到全局极值点,避免陷入局部极值,提高全波形速度反演的精度与速度。

FWI是一个最优化问题,通过对观测波场与理论波场的残差求取极小值来构造目标函数进行迭代更新,最终重构地下介质的参数模型。当前常用的目标函数主要是通过数据残差的l2范数来构建,即:

(1)

我们考虑模型参数m,它对应于速度平方的倒数,是一个实向量m∈N,其中N是空间离散化中的点数。对于每个源和频率,波场和观测数据分别用usv∈N和dobs∈Nr表示,其中,s=1,…,Ns指源,v=1,…,Nv指频率,Nr为接收器的数量。dcal=Pusv,P是将波场投射到接收器位置的算子。

目标函数对m求一阶导数可以得到梯度矩阵,公式为:

(2)

其中,J为Frechet导数矩阵,上标T和H分别表示矩阵的转置和共轭。对目标函数求二阶导数可以得到近似Hessian矩阵,公式为:

Ha=Re[JTJH] 。

(3)

根据高斯牛顿法可以得到模型的更新量为:

(4)

利用梯度下降法[17]最小化目标函数G,模型更新量可以写成:

如果把m表示成N1×N2的速度模型,网格间距为h,用(i,j)代表各个网格点坐标,其中i=1,2,…,N1,j=1,2,…,N2,则在纵向和横向两个方向的梯度分别为:

(6)

全变差定义式如下:

式(7)代表离散模型中每一点的离散梯度的l2范数之和。假设Neumann边界条件使这些差值在边界处为零。我们可以通过定义一个差分算子D来更紧凑地表示‖m‖TV,这样Dm就是离散梯度的一个连接,(Dm)n表示在n索引位置对应的离散梯度向量,n=1,…,N,N=N1×N2,可以定义:

(8)

将全变差约束项应用于常规全波形反演算法中,可以得到基于全变差约束的目标函数,公式为:

对于式(9),如果添加约束m∈[B1,B2]和‖m‖TV≤τ,则式(9)具有等价的优化形式:

mn+Δm∈[B1,B2],‖m‖TV≤τ,

(10)

其中,τ是与λ有关的一个正常数,τ值越小,全变差约束项在目标函数中的权重作用越大。此时反问题可以表示为:

mn+1=mn+Δm

(11a)

Δm∈[B1,B2],‖m‖TV≤τ。

(11b)

根据改进的原始对偶混合梯度(PDHG)方法,求解式(11)的全变差约束问题相当于求解以下问题的鞍点(Δm,p):

(12)

其中,p是拉格朗日乘子,mn+Δm∈[B1,B2],‖·‖∞,2为‖·‖1,2的对偶模,相当于求最大值,而不是l2模的和。

修改后的PDHG方法需要迭代:

mn+Δm∈[B1,B2]

迭代公式可以更明确地写成:

pk+1=pk+δD(mn+Δmk)-

∏‖·‖1,2≤τδ[pk+δD(mn+Δmk)]

(14)

其中,δ和α为迭代步长,满足0<αδ<1/‖DTD‖,∏‖·‖1,2≤τδ(z)表示z在半径为τδ的球上‖·‖1,2范数的正交投影。终止条件为:

(15)

考虑一个2D合成实验,模拟区域网格点大小为200×200,网格宽度h为10 m。图1所示的合成速度模型有一个由较慢的平滑背景包围的恒定的高速区域。将该平滑背景作为反演的初始模型m0。类似于van Leeuwen和Herrmann[18]中的例1,将Ns=34个源置于矩形高速异常体左侧,Nr=81个接收器置于矩形高速异常体右侧。源qsv对应于主频为30 Hz的Ricker小波。

图1 真实速度模型

正则化参数τ有3种不同的选择:①τopt,即慢度平方的总变化量;
②τlarge=1 000τopt, 此时τlarge足够大,全变差约束没有影响;
③τsmall=0.875τopt,比理想的选择稍微小一点。注意,通过使用从式(5)开始的高斯牛顿步骤作为初始模型,式(11)中的凸子问题在τlarge情况下立即收敛。对于所有实验,惩罚参数λ固定为1。从5 Hz数据和6 Hz数据开始,使用计算出的m作为反演6 Hz和7 Hz数据的初始模型,按照由低频到高频的顺序,以此类推。每一个频带都进行50次外部迭代,以确保凸子问题收敛,达到式(15)所示终止条件时,迭代停止。

图2为开展的6次数值实验结果。在无噪声和有噪声的情况下,TV正则化均减少了反演模型中的振荡现象,并使得高速异常体得到更好的恢复。

图2 不同正则化参数τ在有噪声和无噪声情况下反演的速度模型

从噪声数据中反演得到的不连续点的位置要比从无噪声数据中反演的更准确。这是因为噪声在接收器上造成了更大的不连续,从而加强了其他地方TV正则化的效果。τsmall实验证明了这一点,增加TV正则化的权重会导致较大不连续点的分辨率更高。

选取国际标准复杂地质模型开展全波形反演实验。真实速度模型如图3a所示,取自墨西哥湾典型复杂岩体的合成数据——SEAM模型。模拟区域网格点大小为267×384,网格宽度h为20 m。在靠近地面处放置总共72个源和191个接收器。初始模型采用图3b所示平滑模型。图3c及图3d为无噪声数据下的反演结果,图3e及图3f为噪声数据下的反演结果。图3d及图3f为选取τsmall正则化参数下的反演结果,图3c及图3e为选取τlarge正则化参数下的反演结果。由图可知,全变差正则化有助于消除伪影,主要是减少盐下的噪声。

图3 Seam模型的真实速度模型、初始速度模型,以及不同正则化参数τ在有噪声和无噪声情况下反演的速度模型

本文提出了一种在计算上可行的梯度投影算法,该算法可将附加凸约束条件下的全波形反演二次惩罚公式最小化。本文展示了在模型上添加全变差约束和边界约束时,全波形反演中凸子问题的求解过程。数值实验表明,在数据有限或含有噪声的情况下,TV正则化可以通过消除伪影和精确识别不连续点的位置,提高模型反演效果。

本文围绕的梯度投影法,其本质依赖于目标函数的一阶偏微分信息,这使得收敛速度受到一定制约。因此,理论上,我们可以采用二阶方法进行非线性优化,以探索收敛速度是否有提升。

猜你喜欢 变差正则反演 如何深入浅出讲解变差函数的内涵科教导刊(2022年23期)2022-10-15反演对称变换在解决平面几何问题中的应用中等数学(2022年5期)2022-08-29献血后身体会变差?别信!中老年保健(2022年3期)2022-08-24具有逆断面的正则半群上与格林关系有关的同余华南师范大学学报(自然科学版)(2021年5期)2021-11-09反演变换的概念及其几个性质学校教育研究(2018年8期)2018-07-09任意半环上正则元的广义逆上海师范大学学报·自然科学版(2018年3期)2018-05-14sl(n+1)的次正则幂零表示的同态空间华东师范大学学报(自然科学版)(2018年3期)2018-05-14基于ModelVision软件的三维磁异常反演方法地震研究(2017年3期)2017-11-06绿色建筑结构设计指南神州·上旬刊(2017年9期)2017-10-15双次幂变差与价格跳跃的分离金融经济(2014年3期)2015-01-20

推荐访问:反演 梯度 波形

本文来源:http://www.zhangdahai.com/shiyongfanwen/qitafanwen/2023/0605/607423.html

  • 相关内容
  • 热门专题
  • 网站地图- 手机版
  • Copyright @ www.zhangdahai.com 大海范文网 All Rights Reserved 黔ICP备2021006551号
  • 免责声明:大海范文网部分信息来自互联网,并不带表本站观点!若侵害了您的利益,请联系我们,我们将在48小时内删除!