2019, Vol. 62
2. 海洋国家实验室海洋矿产资源评价与探测技术功能实验室, 青岛 266071
2. Laboratory for Marine Mineral Resources, Qingdao National Laboratory for Marine Science and Technology, Qingdao 266071, China
海洋地震勘探中地震波在流-固边界相互耦合的介质(fluid-solid boundary coupled media，简写为FSCM)中传播，对于其数据处理，必须估算FSCM界面处的波传播特性才能准确模拟地震波在FSCM的传播过程.Hou等(2012)和De Basabe和Sen(2015)将地震波在该情况下的数值模拟方法分为两种：(1)单一控制方程法，即在流相和固相中采用相同的地震波波动方程，利用一阶速度-应力方程针对流固耦合介质进行地震波数值模拟.但方程中的变量较多，计算效率相对较低；(2)双方程耦合方法，即分别用声波方程和弹性波方程描述地震波在液相和固相介质中的传播，并利用连续性条件将两者介质相联系.Stephen(1983, 1988)和Sochacki等(1991)利用流固界面的应力、应变连续性条件采用有限差分方法对声学和弹性介质的纯位移方程耦合开展正演模拟研究.Zhang(2012)结合有限差分和有限元方法进行单元格子法FSCM界面地震数值模拟方法的研究以适应海底界面的起伏不规则.Wu和Wu(2017)提出基于高阶有限差分法的一阶速度-应力声波和弹性波方程耦合算法，并模拟验证了地震波在FSCM中传播的稳定性.
全波形反演(Full waveform inversion，简写为FWI)是在最优化理论框架下，通过拟合模拟数据与观测数据，使残差达到最小，从而获取具有高分辨率的地下物性参数(比如纵横波速度、密度、阻抗或各向异性参数等)的反演方法(Virieux and Operto, 2009).Lailly(1983)和Tarantola(1984)最早提出建立最小二乘意义下的目标泛函的FWI框架，利用激发点正传波场和接收点残差反传波场的互相关求取梯度.Pratt和Worthington(1990)将FWI扩展到频率域，并指出多尺度反演能够提高反演的实用性和稳定性.目前FWI已在理论与实际资料中成功应用(Shipp and Singh, 2002; Sears et al., 2008; Sirgue et al., 2010;刘玉柱等, 2014).
多参数FWI是当前研究的热点，但是由于FWI是一个严重的非线性反问题，并且地震数据对不同参数的敏感性在特定的传播角度上有一定的相似性，所以不同参数在反演过程中会相互干涉(Yang et al., 2016a, 2016b).Tarantola(1984)提出基于变密度声波方程的速度和密度的反演策略，Mora(1987)指出近偏移距数据适合阻抗参数反演，而远偏移距数据适合速度参数反演，但两种方法都不能很好地估计密度.Köhn等(2012)指出，在弹性波多参数全波形反演中利用速度-密度参数化反演效果要优于拉梅常数-密度参数化反演.Jeong等(2012)在频率域提出弹性波多参数分步反演策略，首先固定密度反演拉梅常数，然后通过链式法则将结果作为纵横波速度-密度参数化FWI的初始模型，再次反演拉梅常数和密度，以此提高密度反演精度.杨积忠等(2016)将截断高斯牛顿法应用于频率域速度-密度反演，取得不错的效果.对于FSCM下的多参数FWI问题，如果能够通过波形反演精确地估计出P波和S波速度之外的密度，那么对于海洋储层描述和岩性解释将是非常有用的.Wu和Wu(2017)假定密度为常数，提出基于有限差分法的FSCM全波形反演算法，得到良好的纵横波速度反演结果，并对比分析双参数同步反演和分步骤参数反演的效果.
本文提出基于流-固边界耦合变密度介质多参数全波形反演方法.正演部分，在流相介质和固相介质中分别采用二阶声波方程和弹性波位移方程，在流固耦合界面处采用一阶位移-应力弹性波动方程作为流相和固相之间的过渡层，该方法不仅能够提高地震波在流-固边界耦合介质中传播的精度和稳定性，而且极大减少了所需内存，提高计算效率.反演部分，本文提出分步反演策略：第一步，利用海上拖缆接收的地震资料进行FSCM中固相介质P波速度-密度参数化反演，本文基于Plessix(2006)提出的伴随状态法重新推导了基于最小二乘的变密度声波方程反演速度和密度梯度，相比于Tarantola(1984)给出拉梅常数-密度梯度，该方法内存需求更小，同时计算效率更高；第二步，把第一步得到的P波速度结果作为初始速度，进行S波速度参数化反演.这一步反演并未引入密度因素，其一是均匀介质下弹性波位移方程正演模拟效率远远高于非均质方程，其二是密度只对小角度数据敏感，其辐射模式与速度基本一致，所以不考虑密度同样可以得到可信的横波速度.最后，正演方面通过层状模型测试正演的有效性；反演方面对Marmousi-2模型进行数值试验，证明多参数分步全波形反演方法的准确性和适用性.
1 二维流-固边界耦合介质波动方程在各向同性介质情况下，流体中剪切应力为零，一般是由声压参数或位移矢量进行描述，流相中的物理参数一般包括两个：速度和密度.通过一阶速度-压力方程可求得二维声压标量波的波动方程：
|
(1) |
其中，p为声压场，vP为纵波速度，f为震源函数.
在海洋固相介质中地震波是无源的，其传播过程可以由弹性波位移方程描述：
|
(2) |
其中，u, w为位移分量，vS为横波速度.
不难看出，流相和固相控制方程分别对应压力和位移，因此许多学者提出边界连续性条件来模拟地震波的传播(Zhang, 2004；Lee et al., 2009)，
|
(3) |
其中ρ为密度，nx，nz为流-固边界处任意点的法线方向，σxx, σzz为正应力，σxz为切应力.
上述边界耦合算法无法推广到高阶差分格式，随着传播时间增加，界面耦合近似误差积累逐渐增加导致地震波传播不稳定.这里借鉴交错网格的思想，引入一阶位移-应力方程作为过渡层.
在FSCM正演中，本文忽略密度的影响，方便起见令ρ=1，则地震波在弹性介质过渡层中传播满足：
|
(4) |
忽略密度的影响是因为上述方程主要用于S波速度反演中的正演模拟，详见第2节.在流体中方程(4)退化为
|
(5) |
综上所述，通过方程(1)—(5)共同构建基于一阶方程作为过渡层的流-固边界耦合方程(fluid-solid boundary coupled equation，简写为FSCE)，在流相介质由声波控制方程(方程1)模拟，固相介质由弹性波控制方程(方程2)模拟传播，过渡层由一阶位移-应力方程控制(方程4)，包括上层流体和下层固体两部分，两个介质内波场的传播需借助连续性条件(方程5)建立，如图 1所示，具体算法详见附录A.
|
图 1 基于一阶方程过渡层的流-固边界耦合介质正演算法示意图 Fig. 1 Schematic of forward modeling algorithm for medium with fluid-solid coupled boundary based on transition layer of first-order equation |
与传统的交错网格正演模拟相比，本文提出的FSCM正演模拟可以达到与交错网格方法同样的精度阶数.表 1是FSCE与交错网格计算内存的对比，由此可见，FSCE通过减少了方程的变量来降低了计算运行内存，模型纵向深度越大，运行内存占用越小.
|
表 1 流-固边界耦合模型正演模拟内存分析 Table 1 Memory parameters for forward modeling of model with fluid-solid coupled boundary |
Tarantola(1984)给出了正演模拟得到记录ucal与观测记录dobs之间的差异最小的最小二乘法规范意义下的目标函数E：
|
(6) |
其中m表示参数，R表示限定算子.
本文反演采用梯度法，建立迭代格式，
|
(7) |
其中k表示当前迭代轮次，α表示迭代步长，由一维线搜索计算得出，∇mE表示目标函数对参数的梯度.
由于地震数据对不同参数的敏感性在特定的传播角度上有一定的相似性，不同参数在反演过程中会相互干涉，我们把FSCM的多参数FWI分为两个部分：第一步，基于变密度声波方程的P波速度-密度反演；第二步，基于流-固边界耦合方程的S波速度反演.
2.1 基于变密度声波方程的P波速度-密度反演海洋勘探拖缆数据包含了海底固相中丰富的转换波信息，当海底最后通过上覆流体传回拖缆，接收到的都是P波数据，这也对FSCM进行P波速度-密度的双参数反演提供了可行性，这里需要在声波方程引入密度参数作为P波速度-密度反演的控制方程.
首先给出二维非均质下的声波方程：
|
(8) |
De Bartolo等(2012)提出一种基于等效交错网格(Equivalent staggered grid scheme，简写为ESG)的正演算法模拟非均质下的声波方程传播，ESG具有与常规交错网格相同的精度，而且计算效率更高，本文采用ESG法进行正演建模，该方法同时也用于密度梯度的计算，下面会进行阐述.这里将ESG拓展到任意阶差分格式(附录B).
根据Plessix(2006)提出的伴随状态法，梯度构建可以通过正传算子的伴随算子计算得到，由此推出P波速度和密度梯度(附录C)，
|
(9) |
其中，p*表示用伴随算子得到的反传残差波场，∇vPE和∇ρE表示目标函数分别对P波速度和密度的梯度.
Tarantola(1984)给出相同方程下的拉梅常数-密度梯度，通过链式法则得到速度-密度梯度，速度梯度与方程(9)中的一致，密度梯度如下
|
(10) |
对于方程(9)和(10)中密度梯度的时间导数可以省去，这种处理在频率域相当于做了一次低通滤波，根据∇·(αA)=α∇·A+(∇α)A方程(9)可以转换为
|
(11) |
式中的积分项S表示地球表面，若忽略不计即得到和Tarantola一致的密度梯度表达式，因此本文推出的密度梯度更精确一些.而从内存方面看，Tarantola梯度算法需要存储4个波场(正、反传波场p, p*及其导数∇p, ∇p*)信息，而本文给出的梯度算法仅需存储2个波场信息(正、反传波场p, p*)，密度梯度第二项
通过上述方法进行P波速度-密度双参数反演，可以得到准确度高的速度模型，密度模型因为耦合原因结果会差一些，杨积忠等(2016)提出可以先通过双参数反演得到可靠的P波速度，然后利用该结果和原始的密度初始模型重新进行双参数反演，会得到好的结果.
2.2 基于流-固边界耦合方程的S波速度反演为了降低多参数的耦合性，本文在S波速度反演中忽略密度的影响，将第一步得到P波速度作为初始模型，进行S波速度单参数反演.FSCM上覆流体S波速度为零，反演对象限制在海底固相介质，这样就过渡到弹性介质全波形反演.Mora(1987)给出弹性介质的拉梅常数参数化的梯度，
|
(12) |
其中λ和μ表示拉梅常数.通过链式法则得到S波速度梯度，
|
(13) |
通过以上方法最终得到P波、S波及密度参数反演结果，在反演过程中最重要的是控制方程选取，第一步P波速度-密度双参数反演，选取声波变密度方程作为控制方程，利用本文给出的梯度计算得到反演P波速度及密度；第二步S波速度反演过程中，选取FSCE控制方程进行模拟匹配，只更新固相区域S波速度值，最终得到三个参数反演结果.具体流程如图 2所示.
|
图 2 非均匀流-固边界耦合介质多参数全波形反演流程图 Fig. 2 Flow chart of multi-parameter full waveform inversion for heterogeneous medium with fluid-solid coupled boundary |
首先，通过流-固边界耦合的三层模型的数值模拟验证本文正演方法的准确性，我们建立一个300×120网格大小的层状模型，如图 3所示，五角星代表炮点位置(750 m, 5 m)，空间步长统一为5 m，时间采样间隔为0.5 ms，震源选取主频是20 Hz雷克子波，介质为三层介质模型.
|
图 3 流-固边界耦合层状模型 Fig. 3 Layered model with fluid-solid coupled boundary |
图 4展示了分别用声波方程和FSCE正演模拟得到的接收记录，图 4a为纯声波记录，观测到的波有①P0T直达波、②P0TP1RP0T反射波和③P0TP1TP2RP1TP0T反射波三种波，层间也会产生多次波，在记录中振幅很微弱.图 4b为FSCE模拟得到的记录，波形依次为①P0T直达波、②P0TP1RP0T反射波、③P0TP1TP2RP1TP0T反射波、④P0TP1TS2RP1TP0T/P0TS1TP2RP1TP0T反射波和⑤P0TS1TS2RP1TP0T反射波，记录中还存在大量层间反射的转换波，因为能量非常弱而难以显示.上标表示波的类型，T表示透射波，R表示反射波，下标表示波在第几层传播.
|
图 4 流-固边界耦合层状模型地震记录 (a)声波模拟记录；(b) FSCE模拟记录. Fig. 4 Seismic records of layered model with fluid-solid coupled boundary (a) Pure-acoustic record; (b) FSCE record. |
为了更清晰描述FSCE和声波方程模拟的区别，分别抽取FSCE模拟记录和声波模拟记录第50、100和150道进行对比，如图 5所示.红色虚线表示FSCE模拟结果，黑色线条是声波方程模拟的结果.零偏移距时(图 5c)两个记录完全一致；随着偏移距增大(图 5b, 5a)FSCE模拟的结果振幅相对会降低，同时产生一些转换波信息，如图 5a中0.5~0.7 s之间，而采用声波方程模拟则得不到这些结果.
|
图 5 流-固边界耦合层状模型地震记录 (a)第50道；(b)第100道；(c)第150道；(红色虚线是FSCE模拟结果，黑色线条表示声波方程模拟结果). Fig. 5 Seismic records of layered model with fluid-solid coupled boundary (a) 50th trace; (b) 100th trace; (c) 150th trace; Red lines are the FSCE, black lines are the pure-acoustic. |
表 2为本文提出的过渡层法和传统交错网格差分方法运行内存和计算时间对比结果.通过对比两种方法的内存使用情况，可见本文方法降低了将近一半正演模拟的运行内存，通过表 1分析可知，模型纵深增大，内存优势更加明显.另一方面，过渡层法只在流-固耦合处采用交错网格剖分，在液相和固相中均使用二阶方程，通过对比两种方法对流固单层模型的模拟时间可见，本文方法在计算效率上有所提高.
|
|
表 2 流-固边界耦合层状模型正演模拟内存及运行时间 Table 2 Memory and operation time of forward modeling on layered model with fluid-solid coupled boundary |
为测试本文反演方法的准确性以及在地表观测系统下反演复杂模型能力，本节采用Marmousi-2模型进行测试.模型大小为2400 m×1360 m，其中水深为120 m，网格间距为8 m，震源采用主频为20 Hz的雷克子波，时间采样间隔为0.4 ms，震源激发位置位于水中坐标(1200 m, 8 m)的位置，检波器布于水面.
图 6展示了真实的纵横波速度(图 6a, 6b)以及密度模型(图 6c)，通过高斯平滑后得到初始速度(图 6d, 6e)和密度模型(图 6f).在三参数反演中，首先利用观测数据进行纵波速度和密度反演，本文采用同时反演的策略，迭代次数为80次，用反演得到的纵波速度作为第二步横波速度反演中的初始模型，迭代次数为50次，反演结果如图 7所示.图 8显示了反演结果的纵向剖面.通过对比可以看出，纵波和横波速度反演结果具有较高的分辨率与反演精度，而密度反演结果更像是偏移剖面，反演效果与真实模型一致性稍差，尤其是水平方向1200~1600 m，深度1000 m左右的地质构造反演较差，其他区域结果比较理想.这是因为密度辐射模式能量集中于小角度，只需小偏移距的反射波即可快速反演，其收敛速度比纵波速度要快一些，这种不同参数收敛速率的差异性，导致密度反演结果变化更剧烈，而速度反演结果更平滑.
|
图 6 Marmousi-2模型 (a)纵波速度真实模型；(b)横波速度真实模型；(c)密度真实模型；(d)纵波速度初始模型；(e)横波速度初始模型；(f)密度初始模型. Fig. 6 Marmousi-2 models (a) Real P-wave velocity model; (b) Real S-wave velocity model; (c) Real density model; (d) Initial P-wave velocity model; (e) Initial S-wave velocity model; (f) Initial density model. |
|
图 7 FSCM多参数全波形反演结果 Fig. 7 Results of FSCM multi-parameter full waveform inversion |
|
图 8 多参数全波形反演结果纵向剖面 (a)纵波速度，(b)横波速度，(c)密度，图中的5条线依次位于水平方向250 m、500 m、750 m、1000 m、1250 m处(黑色线条是模型真实值，黑色虚线是反演的初始值，红色实线是最终反演结果). Fig. 8 Vertical profiles at 250 m, 500 m, 750 m, 1000 m and 1250 m in horizontal direction of the models shown in Fig. 6 and Fig. 7 (a) P-wave velocity; (b) S-wave velocity; (c) density; Dotted lines are initial models, black lines are true models, red lines are the final inversion. |
下面通过对比弹性波多参数反演验证本文方法的优势，需注意FSCM的特征是上覆流体，即横波零速层，当速度为零会造成数值模拟不稳定.为满足稳定性条件，这里将Marmousi-2模型上覆流体的横波速度设置为1000 m·s-1，建立的弹性层状介质模型即可由方程(2)实现模拟.图 9显示的是FSCM与纯弹性介质下分别使用FSCE与弹性波方程模拟的波场，事实上，FSCE模拟海洋勘探(图 3)接收的一定是纯P波资料，其震源激发后在海底的一次反射(图 9a)，来自海底固相层界面的反射波再次穿过固-流耦合边界到达检波点(图 9c)都不会产生转换S波，虽然在弹性波正演模拟可以加载纵波源模拟海上震源激发，但是遇到海底反射(图 9b)以及上行波穿过海底界面(图 9d)均会产生转换S波.再结合图 4中与标量波模拟结果分析对比，能够真实还原海洋环境是FSCE的优势所在.
|
图 9 FSCM和弹性介质波场对比 (a) t=0.3 s时刻FSCM波场; (b) t=0.3 s时刻弹性介质波场; (c) t=0.5 s时刻FSCM波场; (d) t=0.5 s时刻弹性波场. Fig. 9 Wavefields of FSCM and elastic medium (a) 0.3 s, FSCM wavefield; (b) 0.3 s, elastic medium wavefield; (c) 0.5 s, FSCM wavefield; (d) 0.5 s, elastic medium wavefield. |
需要指出，图 9a和图 9c的波场是分为流相和固相两部分，流相中的波场由压力p控制，而固相中由位移场u，w控制(图中展示压力p-水平分量位移u波场)，二者量纲不同，所以对位移乘以系数以得到可视化波场.统一起见，图 9b和图 9d展示为正应力τzz-水平分量位移u波场.
下面通过与弹性波多参数同步反演对比，验证非均匀流-固边界耦合介质多参数全波形反演方法.与弹性层状模型类似，我们需要先将Marmousi-2模型上覆流体S波速度设为800 m·s-1，然后进行三参数同步反演.图 10为反演结果，显然vP和vS反演目标层界面偏离真实位置(黑色圆圈部分)，而ρ的反演结果已经完全偏离了真实模型，一是因为反演过程中不同参数相互耦合，二是多参数同步反演非线性太强.抽取反演结果的纵向剖面，与本文提出的反演方法比较(图 11).通过对比可以看出，弹性波三参数同步反演得到的vP和vS只存在更新界面，而速度层析很差的问题(黑色圆圈部分)，也是由目标层反演偏离所造成的.反观流-固边界耦合介质多参数反演结果要优于弹性波反演.当然弹性波多参数分步反演策略(Jeong et al., 2012)已经取得很好效果，本文采用弹性波同步反演原因在于非均匀流-固边界耦合介质具有上覆横波零速层的特点，数值模拟中不能用弹性波方程直接进行反演，同时其观测记录亦不包含S波.因此本文提出的方法与弹性波多参数同步反演本不具备对比性，但是从反演策略角度分析，通过两者对比也能证明非均匀流-固边界耦合介质多参数全波形反演方法相对于传统方法效果更好.
|
图 10 弹性波三参数全波形反演结果 Fig. 10 Results of elastic wave multi-parameter full waveform inversion |
|
图 11 多参数全波形反演结果纵向剖面 (a)纵波速度; (b)横波速度; (c)密度，图中的5条线依次位于水平方向250 m、500 m、750 m、1000 m、1250 m处(黑色线条是真实模型，黑色虚线是反演的初始值，红色实线是FSCM最终反演结果，蓝色实线是弹性介质最终反演结果). Fig. 11 Vertical profiles at 250 m, 500 m, 750 m, 1000 m and 1250 m in horizontal direction of the models shown in Fig. 10 (a) P-wave velocity; (b) S-wave velocity; (c) density; Dotted lines are the initial model, black lines are the true model, red lines are the final inversion of FSCM, blue lines are the final inversion of elastic medium. |
本文针对流-固边界耦合介质提出了一种高效、稳定的多参数全波形反演方法，并通过模型测试得出以下结论：
(1) 本文采用一阶位移-应力程作为过渡层连接声波方程和弹性波位移方程来模拟海洋勘探环境，该方法可以推广到高阶有限差分建模，解决的采用边界条件耦合导致的不稳定性，同时相比于交错网格建模，大大降低了计算成本和内存.
(2) 反演策略对多参数全波形反演至关重要，本文将其分为先反演纵波速度和密度，再反演横波速度的两个步骤来完成.文中利用伴随状态法得到新的密度梯度并采用等效交错网格法求解，以计算代替存储，减少了储存量.通过Marmousi-2模型的数值试验测试说明反演方法的准确性和适用性.
致谢感谢国家科技重大专项(2016ZX05024-001-008)的资助，特别感谢评审专家为本文完善提出的宝贵意见.
附录A(FSCM边界耦合算法)借助交错网格可以稳定数值模拟地震波在流固界面传播的特点，在流固界面处采用一阶位移-应力波动方程差分模拟，与上覆流体中的声压标量波方程和下伏固相中的纯位移弹性波方程的耦合.
如附图 1所示，假设流-固界面为水平界面，N_fluid为流-固耦合界面的深度，N_fluid－5到N_fluid+5是交错网格区域.首先利用方程(A1, A2)计算n-1时刻N_fluid－15到N_fluid+10之间的应力，然后用一阶位移-应力方程更新n时刻N_fluid－10到N_fluid+5之间的位移，然后利用方程(A2)计算n时刻N_fluid－5到N_fluid的应力，最后通过声波方程和弹性波方程分别更新n时刻1到N_fluid－6之间的压力以及n时刻N_fluid+6到Nz之间的位移.
|
(A1) |
|
(A2) |
|
附图 1 使用一阶位移-应力方程过渡FSCM耦合算法示意图(浅色区域是流体，深色区域是固体) Appendix Fig. 1 FSCM coupled scheme with the transition of first-order displacement-stress equation (light-colored area is fluid and dark-colored area is solid) |
为实现无限远模拟区域的假设，采用完全匹配层吸收边界条件，即在液相和固相中采用二阶常规PML吸收边界条件，在过渡层中采用一阶常规PML吸收边界条件.
附录B(等效交错网格任意阶差分格式)传统的交错网格差分格式是在空间和时间上应力和速度交错更新，每次都需要保存所有参数值，而等效交错网格在利用其周围应力或应变信息时，又可以采用其本身的差分格式来代替，从而实现仅含有纯压力(应力)或纯位移交错差分思想，如附图 2所示.
|
附图 2 交错网格与等效交错网格对比 左图交错网格剖分；右图准规则网格剖分. Appendix Fig. 2 Staggered grid and equivalent staggered grid scheme The left panel shows the staggered grid, the right panel shows the equivalent staggered grid. |
De Bartolo等(2012)在解决非均质声波方程正演问题中给出ESG四阶差分格式，这里将其推广到任意阶精度，下面给出空间任意阶差分格式：
|
(B1) |
其中n表示当前时刻，ai表示一阶导数的差分系数，b在非均质声波方程中表示密度的倒数，而在本文中应用于梯度计算，相当密度平方的倒数，dx表示网格间距，N表示差分阶数，对于z方向与之类似.本文使用的是空间10阶差分精度，接下来给出差分格式(B2)以及差分系数(附表 1).
|
(B2) |
|
|
附表 1 任意阶差分系数 Appendix Table 1 Arbitrary order difference coefficient |
在非线性反演过程中，计算量最大的部分为梯度的计算，为此很多学者做了大量的研究.Plessix(2006)引入共轭状态法求解梯度，并在范函分析的意义下进行了详尽的数学定义.
时间域声波变密度方程表达如下:
|
(C1) |
建立状态方程，
|
(C2) |
对目标函数引入参数扰动，
|
(C3) |
其中λ为残差的伴随波场，〈·〉表示内积，利用文献中(Plessix, 2006)的增广函数法求出伴随方程，
|
(C4) |
其中p*(t)=λ(T－t)残差反传时间是从最大时间到零时刻.将方程(C4)代入目标函数的扰动方程(C3)，
|
(C5) |
由此计算得到速度和密度的梯度表达式，
|
(C6) |
在式(C6)中速度和密度的梯度可以做一步近似，即略去对时间的二阶导数项，相当于在频率域做了一次低通滤波，在实际求解时，能够起到避免迭代出现局部极小值的效果.
De Basabe J D, Sen M K. 2015. A comparison of finite-difference and spectral-element methods for elastic wave propagation in media with a fluid-solid interface. Geophysical Journal International, 200(1): 278-298. DOI:10.1093/gji/ggu389 |
Di Bartolo L, Dors C, Mansur W J. 2012. A new family of finite-difference schemes to solve the heterogeneous acoustic wave equation. Geophysics, 77(5): T187-T199. DOI:10.1190/geo2011-0345.1 |
Hou G N, Wang J, Layton A. 2012. Numerical methods for fluid-structure interaction-a review. Communications in Computational Physics, 12(2): 337-377. DOI:10.4208/cicp.291210.290411s |
Jeong W, Lee H Y, Min D J. 2012. Full waveform inversion strategy for density in the frequency domain. Geophysical Journal International, 188(3): 1221-1242. DOI:10.1111/j.1365-246X.2011.05314.x |
Köhn D, De Nil D, Kurzmann A, et al. 2012. On the influence of model parametrization in elastic full waveform tomography. Geophysical Journal International, 191(1): 325-345. DOI:10.1111/j.1365-246X.2012.05633.x |
Lailly P. 1983. The seismic inverse problem as a sequence of before stack migrations.//Conference on Inverse Scattering: Theory and Application.Philadelphia, PA: Society for Industrial and Applied Mathematics, 206-220.
|
Lee H Y, Lim S C, Min D J, et al. 2009. 2d time-domain acoustic-elastic coupled modeling:a cell-based finite-difference method. Geosciences Journal, 13(4): 407-414. DOI:10.1007/s12303-009-0037-x |
Liu Y Z, Xie C, Yang J Z. 2014. Gaussian beam first-arrival waveform inversion based on Bornwavepath. Chinese Journal of Geophysics (in Chinese), 57(9): 2900-2909. DOI:10.6038/cjg20140915 |
Mora P. 1987. Nonlinear two-dimensional elastic inversion of multioffset seismic data. Geophysics, 52(9): 1211-1228. DOI:10.1190/1.1442384 |
Plessix R E. 2006. A review of the adjoint-state method for computing the gradient of a functional with geophysical applications. Geophysical Journal International, 167(2): 495-503. DOI:10.1111/j.1365-246X.2006.02978.x |
Pratt R G, Worthington M H. 1990. Inverse theory applied to multi-source cross-hole tomography. Part 1:Acoustic wave-equation method. Geophysical Prospecting, 38(3): 287-310. DOI:10.1111/j.1365-2478.1990.tb01846.x |
Sears T J, Singh S C, Barton P J. 2008. Elastic full waveform inversion of multi-component OBC seismic data. Geophysical Prospecting, 56(6): 843-862. DOI:10.1111/j.1365-2478.2008.00692.x |
Shipp R M, Singh S C. 2002. Two-dimensional full wavefield inversion of wide-aperture marine seismic streamer data. Geophysical Journal International, 151(2): 325-344. DOI:10.1046/j.1365-246X.2002.01645.x |
Sirgue L, Barkved O I, Dellinger J, et al. 2010. Thematic set:fullwaveform inversion:the next leap forward in imaging at Valhall. First Break, 28(4): 65-70. |
Sochacki J S, George J H, Ewing R E, et al. 1991. Interface conditions for acoustic and elastic wave propagation. Geophysics, 56(2): 168-181. DOI:10.1190/1.1443029 |
Stephen R A. 1983. A comparison of finite difference and reflectivity seismograms for marine models. Geophysical Journal International, 72(1): 39-57. DOI:10.1111/j.1365-246X.1983.tb02803.x |
Stephen R A. 1988. A review of finite difference methods for seismo-acoustics problems at the seafloor. Reviews of Geophysics, 26(3): 445-458. DOI:10.1029/RG026i003p00445 |
Tarantola A. 1984. Inversion of seismic reflection data in the acoustic approximation. Geophysics, 49(8): 1259-1266. DOI:10.1190/1.1441754 |
Virieux J, Operto S. 2009. An overview of full-waveform inversion in exploration geophysics. Geophysics, 74(6): WCC1-WCC26. DOI:10.1190/1.3238367 |
Wu J H, Wu G. 2017. Full waveform inversion in Acoustic-elastic coupled media based on the finite-difference method.//87th Ann. Internat Mtg., Soc. Expi. Geophys.. Expanded Abstracts.
|
Yang J Z, Liu Y Z, Dong L G. 2016. Multi-parameter full waveform inversion for velocity and density based on Born sensitivity kernels. Chinese Journal of Geophysics (in Chinese), 59(3): 1082-1094. DOI:10.6038/cjg20160328 |
Yang J Z, Liu Y Z, Dong L G. 2016a. Least-squares reverse time migration in the presence of density variations. Geophysics, 81(6): S497-S509. DOI:10.1190/geo2016-0075.1 |
Yang J Z, Liu Y Z, Dong L G. 2016b. Simultaneous estimation of velocity and density in acoustic multiparameter full-waveform inversion using an improved scattering-integral approach. Geophysics, 81(6): R399-R415. DOI:10.1190/geo2015-0707.1 |
Zhang J F. 2004. Wave propagation across fluid-solid interfaces:a grid method approach. Geophysical Journal International, 159(1): 240-252. DOI:10.1111/j.1365-246X.2004.02372.x |
刘玉柱, 谢春, 杨积忠. 2014. 基于Born波路径的高斯束初至波波形反演. 地球物理学报, 57(9): 2900-2909. DOI:10.6038/cjg20140915 |
杨积忠, 刘玉柱, 董良国. 2016. 基于Born敏感核函数的速度、密度双参数全波形反演. 地球物理学报, 59(3): 1082-1094. DOI:10.6038/cjg20160328 |