《层析成像中的波形反演》由会员分享，可在线阅读，更多相关《层析成像中的波形反演（6页珍藏版）》请在金锄头文库上搜索。
1、第 18 卷 第 4 期 1996 年 11 月 物探化探计算技术 Vol. 18 No. 4 Nov. 1996 COM PUT INGT ECHNIQUESFORGEOPHYSICALANDGEOCHEMICALEXPLORAT ION层析成像中的波形反演?王大川? ?严忠琼 李显贵( 成都理工学院)【 摘 要】 在波动方程有限差分数值模拟基础上, 用波形记录残差平方并建立目标函数。对目标函数进行梯度优化, 我们得出了介质速度的迭代公式。通过计算目标函数梯度我们可以对介质进行速度校正, 从而获得反演解。经数值试验证明这是一种很有效的反演方法。【 关键词】 层析成像 波形反演 梯度法0 引言
2、地震层析成像由于能提供较高的分辨率而在地球物理界受到极大关注, 成为当今地球物理学最前沿的研究课题之一。 由于射线理论存在明显的局限性, 所以今后波动方程层析成像将成为地震层析成像的主流。1 有限差分正演模拟解波动方程的数值方法主要有有限元素法、 有限差分法、 Fourier 变换法和积分法等。有限差分法通过对波场函数的二阶偏导作差分逼近推导出波动方程的差分格式。有限差分正演计算容易在计算机上实现, 且具有较快的速度和较高的精度, 所以在波形反演中广泛使用。二维声波方程可写为:1 V2?2U?t2- ?2U= f( 1)其中: U= U( x , z, t) 为波场函数; V = V ( x,
3、 z ) 为空间速度分布; ?2=?2?t2?2?z2为空间二维Laplace 算子; f = f ( xs, zs, t) 为震源函数。应用差分理论将波场函数的递推公式表示为:U( x, z , t+ ?t) = 2U( x , z, t) - U( x, z , t- ?t) +V2?t2?x2Tx +V2?z2?tzT z+ V2?t2f ( 2)其中T x= U( x + ?x, z, t) - 2U( x, z, t) + U( x - ?x, z, t)T z= U( x , z+ ?z, t) - 2U( x, z, t) + U( x , z- ?z, t)?收稿日期: 199
4、6- 05- 02; 改稿日期: 1996- 08- 16 地址: 成都市二仙桥东三路 1 号( 邮政编码: 610059)为地矿 85 重点研究项目“ 勘探地球物理层析成果技术研究”对于两个较典型的模型作井间正演模拟计算, 取其井间距为 30m, 计算井段为 60m, 划分为 3161 个网格, 即网格边长为 1m。 时间采样间隔为 0. 1ms, 记录长度为 100ms。 观测系统布置如下图所示。 3 个炮点等距排列, 61 个接收点在另一井等距排列。 子波函数采用频率为 150Hz 的雷克子波, 每个模型均获得三个共炮点记录。模型2模型1观测系统- - 炮点- - 接收点v3= 4000
5、m/ sv2= 3000m/ sv1= 2000m/ sv1= 3000m/ sv3= 3000m/ sv2= 2000m/ s40202040606000303030060模型 1 的三个共炮点记录如图 1 所示, 记录中可以明显地观察到透射波及界面反射波。图 1 三个共炮点记录图 2 的模拟记录中, 界面的反射波也很明显, 在第二个共炮点记录中, 可以清楚地看到导波沿低速层传播的轨迹。图 2 模拟记录2994 期王大川等: 层析成像中的波形反演由上述两模型的正演模拟计算结果可以看出: 有限差分正演能很好地模拟地震波在 介质中的传播。2 波动方程梯度法反演反演方法可大体划分为线性和非线性两大
6、类。线性反演一般对地球物理问题作了一定程度的线性简化, 使之退化为求解线性方程组的问题。 梯度法是一种常用的函数优化方法, 属于非线性方法的范畴, 它不需要求解矩阵方程。 下面就讨论波动方程梯度反演方法。2. 1 目标函数与迭代公式时域波形反演的目标函数写为记录残差平方和的形式, 即:? ( V) =1 2srT0( U)2dt( 3)其中: U= Uobv- Ucal; Uobv为观测记录值; Ucsl为计算记录值; T 为处理记录长度; V = ( V1,V2, , Vn) 为介质速度( 将介质区域离散为 n 个网格节点) 。将 ? 进行 T aylor 展开并略去二阶及二阶以上导数项:?
7、 ( V+ ?V ) = ?( V) +i? ?Vi?V?Vi( 4)反演的过程即是极小目标函 数的过程, 所以要求目标函数值逐次下降, 即要求:?( V ) - ?( V+ ?V) 0 由式( 4) 可得:?( V) - ?( V + ?V ) = - (? ?V1?Vi+? ?V2?V2+ +? ?Vn?Vn)V= - (? ?V1,? ?V2, ,? ? Vn) ( ?V1, ?V2, , ?Vn)T= - g( V ) ?VT( 5)其中 g( V ) = (? ? V1,? ?V2, ,? ?Vn) 是目标函数在 V 的一阶导数的各分量, 称为目标函数的 梯度或斜量。显然, 要使式(
8、 5) 有最大值, 则应使点积:- g( V ) ?V = - gV?V cos!( 6)有最大值。 其中 !为 g 与 ?V 之间的夹角。 显然当 ! = 时, 即进行寻查极小值方向 ?V 与 梯度方向 g 相反时, 才能保证最快速下降, 沿目标函数负梯度方向进行单向寻查, 求得下降最低的点。故迭代公式写为:V = V - # g( 7)其中 #为最优步长。如此反复迭代, 直到目标函数 ?( V ) 达到预定精度为止。2. 2 梯度的求取 由式( 3) 对 V 求导得梯度表达式:g=? ?V1=sr?( U) ?V U=sr?Ucal ?V U( 8)其中?Ucal ?V为 Frecher
9、导数。对声波方程, 相应的格林函数满足:1 V2?2G? t2- ?2G= ( 9)式中: G= G( x, z , t, x , z , t ) 为 Green 函数; = ( x , z , t ) 为时空 函数。Green 函数可理解为位于( x , z ) 的 t 时刻的集中体力 ( x , z , t ) 在速度介质空间上产生300物探化探计算技术18 卷的波场函数。 根据Green 公式, 可将声波方程的解表示为 Green 函数积分形式:U( x , z, t) =t0dt?sdx dz G( x , z, t, x , z , t ) f ( xs, zs, t )根据褶积的定
10、义及 Green 函数的互易性, 可将上式写为:U( x , z, t) =?sdx dz G( x, z, t, x , z , 0) * f ( xs, zs, t)( 10)假设速度扰动 V 引起的散射场为 U, 则( V+ V) 和( U+ U) 同样满足声波方程:1 ( V + V )2?2( U+ U)?t2- 2( U+ U) = f( 11)采用下列近似:1 ( V + V )2=1 V2-2 V V3代入式( 11) 与式( 1) 相减, 整理后得:1 V2?2( U)?t2- 2( U) =?2( U+ U)?t22 VV3( 12)右边采用弱散射近似, 即令 U+ U!
11、U, 则式( 12) 写为:1 V2?2( U)?t2- 2( U) =?2U?t22 VV3( 13)式( 13) 的解为: U=?sdx dz G* (?2U?t22 VV3)( 14)由于只考虑一个点上的速度扰动, 所以上式写为:U= G* (?2U?t22 VV3) =2 V V2G*?2U?t2两边同除以 V , 则可得出Frecher 导数的表达式:?Ucal ?t=2 V3G*?2U?t2所以梯度的表达式写为:g=sr2 V3G*?2U?t2( 15)观测系统图3炮点 接收点010103 反演数值计算实例正演模拟计算采用有限差分法, 震源波子采用150Hz 的雷克子波。成象区域划
12、分为 1111 个网 格, 每个网格边长为 1m, 时间采用间隔为 0. 1ms, 时间采集长度为 10ms。观测系统布置如图 3 所示, 在左井中等距布置着 11 个炮点, 11 个接收点在另一井中等距放置, 共有 121 个记录。 模型 1: 理论模型如图 4所示, 由速度为 3000m/ s的均匀模型经15 次迭代后得到图 5 的结果。 可以看出,高速异常基本得到恢复, 只是数值上略有差异。模型 2: 理论模型如图 6 所示, 初始迭代选用速度为 3000m/ s 的均匀模型, 20 次迭代3014 期王大川等: 层析成像中的波形反演后得到图 7 所示结果。 在炮点所在井附近, 速度恢复
13、好, 而在远离炮点的地方, 速度恢复不 是太理想, 所以计算梯度时要将到炮点的距离考虑进去。图 6 理论模型图 7 20 次迭代后结果图模型 3: 这是一个复杂的断层模型, 如图 8 所示。初始迭代也采用速度为 3000m/ s 的均匀模型, 经 20 次迭代后的结果如图9 所示。 可见断层的形状基本恢复, 只是速度相差较 大,图 8 断层模型图 9 20 次迭代后结果图302物探化探计算技术18 卷这是由于理论模型的速度变化值太大( 从 2000m/ s 到 5000m/ s) 造成的, 所以波动方程的 反演要求介质的速度扰动不能过大。4 结束语虽然数值试验的结果预示着波动方程层析成像具有较
14、大的发展潜力, 但应该认识到对波动层析的研究仍处在理论阶段。 弱散射近似降低了方法的实用性, 现场资料的复杂性以及提取某些参数( 像震源子波) 的艰巨性都将严重地阻碍波动层析的发展应用。参考文献1 Odile Gauthier.Jean Virieux,and Albert T arantola.T wo-dimension nonlinear inversion of seismicwaveforms:Numerical results. Geophysics. 1996; 51: 1387- 14032 Woodward. M J. Wave- equation tomography. G
15、eophysics. 1992; 57: 15- 263 Alford. R M . Kelly. and Boore. D M . Accurancy of finite difference modeling of the acoustic wave- equation.Geophysics. 1974; 39: 834- 8434 Albert C.Beynolds.Boumdary condituins for the numerical solution of wave propagation problems.Geophysics. 19; 43( 6) : 1099- 11105
16、 杨文采. 地球物理反演和地震层析成像. 北京: 地质出版社.19896 朱介寿等. 地震学中的计算方法. 北京: 地震出版社. 19887 尹峰等. 最大熵时域逆散射层析成像的研究. 地球物理学报. 1994; 37: 673- 6818 ( 美) 安艺敬一.P G. 理查兹. 定量地震学理论与方法. 北京: 地震出版社. 1986WAVE- FORM INVERSION FOR TOMOGRAPHYWang Dachuan Yan Zhongqiong Li Xiangui (Chengdu Institute of Technology)Abstract The objective function can be built based on the finite difference modeling
