地球物理学报  2020, Vol. 63 Issue (8): 3078-3090   PDF    
地震波模拟中多轴复频移近似完全匹配层的吸收效果及稳定性
罗玉钦, 刘财     
吉林大学地球探测科学与技术学院, 长春 130026
摘要:在进行地震波模拟计算的过程中用有限的计算区域模拟地下无限空间,需要进行边界截断.为了在边界处不产生虚假反射影响模拟结果,需要引入吸收边界条件.本文采用的近似完全匹配层是一种新型非分裂完全匹配层,计算效率较高.同时相比于其他非分裂完全匹配层,其还具有不改变方程的形式、易于实现等优势.但是当入射波角度较大,边界吸收效果变弱,且残留在边界中的能量使近似完全匹配层变得极其不稳定.多轴复频移近似完全匹配层的提出就是为了改善对大角度入射波的吸收并且提高边界的稳定性.通过实验模拟和矩阵特征值灵敏度来研究多轴复频移近似完全匹配层的吸收效果及稳定性.结果表明该方法不仅能够吸收掠入波,而且对常规入射波的吸收也得到提升,同时拥有更好的稳定性.
关键词: 完全匹配层      地震波模拟      多轴      复频移     
On the stability and absorption effect of the multiaxial complex frequency shifted nearly perfectly matched layers method for seismic wave propagation
LUO YuQin, LIU Cai     
College of Geo-exploration Science and Technology, Jilin University, Changchun 130026, China
Abstract: In the process of seismic wave simulation calculation, a limited calculation area is used to simulate the underground infinite space, and boundary truncation is required. In order to prevent false reflections at the boundary which will affect the simulation results, it is necessary to introduce absorption boundary conditions. The nearly perfectly matched layer (NPML) used in this paper is a new type of non-split perfect matched layer with high computational efficiency. At the same time, it does not change the form of the equation and is easier to implement than other non-split perfectly matched layers. However, when the incident wave angle is large, the boundary absorption effect becomes weak, and the energy remaining in the boundary makes the nearly perfect matched layer extremely unstable. The multiaxial complex frequency shifted-NPML (MCFS-NPML) is proposed to improve the absorption of large angle incident waves and the stability of the boundary. The absorption effect and stability of MCFS-NPML were studied by experimental simulation and matrix eigenvalue sensitivity. The results show that the method can not only absorb the grazing wave, but also improve the absorption of the normal incident wave, and has better stability.
Keywords: Perfectly matched layers    Simulation of seismic wave    Multiaxial    Complex frequency shifted    
0 引言

Berenger(1994)提出将完全匹配层作为电磁波传播问题中的边界条件.完全匹配层围绕着物理区域,将来自物理区域中的任何方向和频率的波吸收衰减.Chew和Weedon(1994), Collino和Monk(1998)将其解释为复坐标伸展变换的结果.完全匹配层不仅在麦克斯韦方程中运用,很快扩展到亥姆霍兹方程(Teixeira and Chew, 1998; Harari et al., 2000)、线性欧拉方程(Hu, 1996;Nataf, 2005)、弹性波动方程(Chew and Liu, 1996; Collino and Tsogka, 2001)等情况之中.

随着研究的深入,发现该方法在入射波角度很大时,衰减系数会变得很小,因此传统的完全匹配层无法吸收掠入波,在极低频时会产生奇异值.在模拟过程中发现完全匹配层是不稳定的.Kuzuoglu和Mittra(1996)对复坐标伸展变换形式进行修改,原来为1的项被尺度因子替换,在分母中引入频移因子.频移因子可以防止在极低频时产生奇异值,尺度因子则在对大角度入射波的吸收中起作用.随后Roden和Gedney(2000)提出了基于该方法下的完全匹配层——复频移卷积完全匹配层(C-PML),Komatitsch和Martin(2007)将C-PML推广到弹性介质中,同时指出该边界的不稳定性.Festa等(2005)指出PML边界对面波的吸收能力较弱且会产生不稳定现象.Meza-Fajardo和Papageorgiou(2008)认为传统的PML在各向同性和各向异性中都不满足严格的渐进稳定性,并引入多轴技术,即在多个正交方向上同时引入衰减因子.近年来国内的学者对完全匹配层做了许多研究(王守东, 2003; 刘洋, 2014; 廉西猛等, 2015):将完全匹配层扩展到探地雷达模拟(冯德山等, 2016; 冯德山和王珣, 2017)、全波形反演(景旺等, 2018)、伪谱法模拟(黄建平等, 2016)、瞬变电磁法模拟(李展辉和黄清华, 2014)、大地电磁模拟(薛帅等, 2017)等方面的应用.田坤等(2014)将多轴技术应用于C-PML中;谢志南和章旭斌(2017)提出采用复坐标延伸技术变换弱形式波动方程构建完全匹配层;秦臻(2009)对PML进行改进,熊章强和毛承英(2011)在秦臻的基础上引入复频移变换;张显文等(2009)对复频移递归积分PML进行研究;陈可洋(2010)提出三角函数型衰减函数,罗玉钦和刘财(2018)设计的衰减函数能够削弱离散带来的虚假反射,对衰减因子进行了系统的研究.

与其余非分裂完全匹配层不同,近似完全匹配层NPML(nearly PML)是直接对波场变换求取添加边界之后的方程组,方程组的形式没有改变,不需要卷积处理.由Cummer(2003)提出,随后Berenger(2004)对NPML的吸收效果进行研究,很快被推广到声波介质、弹性介质以及双相介质中(Hu et al., 2007; Chen, 2011, 2012).该方法同样存在着吸收与稳定性方面的问题.因此我们提出多轴复频移近似完全匹配层(MCFS-NPML),保留近似完全匹配层的优势,同时改善边界对掠入波的吸收效果并通过引入频移因子与稳定性因子提高边界稳定性.引入的频移因子、尺度因子和稳定性因子是MCFS-NPML在改善NPML吸收效果和稳定性的关键,本文通过求解方程组,分析矩阵特征值以及模拟结果对各个因子进行系统研究,也证明了MCFS-NPML拥有更好的吸收效果和稳定性.

1 理论研究 1.1 近似完全匹配层

假设外力为0,弹性动力学问题可以通过柯西方程和胡克定律来描述:

(1)

u是位移场,τ是应力张量,E是应变张量且E=1/2[u+(▽u)T],C是弹性系数张量,ρ是密度.速度是位移对时间求取一阶导,用速度变量将式(1)化为一阶速度应力方程:

(2)

为了推导出近似完全匹配层下的新方程组,将式(2)变换到频率域中:

(3)

式中:ω是圆频率;V为频率域弹性波场速度分量;T为频率域弹性波场应力分量.采取的复坐标伸展变换(CCS)的形式为

(4)

(5)

将式(5)代入式(3),并直接作用在速度分量和应力分量上,得到复坐标伸展变换后的方程组:

(6)

联立式(5)式(6)有(以Txx为例)

(7)

(8)

同理,对其他的变量采取相同的方法并变换到时间域中:

(9)

其中变换之后的波场可以通过(10)式求得

(10)

通过式(9)我们可以看到主控方程的形式并没有发生变化.这里需要引入8个辅助变量来替换原主控方程中的变量,通过式(10)中的微分方程来求取.在计算区域中衰减因子值为0,因此在编程的时候可以只在PML区域中进行波场变换来减少存储空间.

1.2 复频移近似完全匹配层

通过式(4)对波场进行复拉伸变换得到近似完全匹配层,将复拉伸变换公式做如下改动:

(11)

对比式(11)和式(5),等式右侧的1变为与x相关的函数,称之为尺度因子;分母中引入频移因子.采用公式(11)对式(2)做变换,过程如下:

(12)

(13)

转换到时间域:

(14)

对其余变量采取相同的处理方法得到变换公式:

(15)

用公式(15)替换公式(10)得到引入频移因子与尺度因子后的变换形式.

1.3 多轴复频移近似完全匹配层

完全匹配层的分布如图 1所示,在区域1,αz=αz(z)αx=0.该区域衰减因子在主方向上呈指数增长,但在与主方向垂直的方向上衰减因子值为0,相当于只有一个剖面,而在区域3,右侧与顶部的衰减剖面重合,在角区间存在双衰减剖面,因此我们引入稳定性因子,增加新剖面的时候只对区间1区间2处理即可.

图 1 衰减因子PML区域示意图 Fig. 1 The schematic of attenuation coefficients in the PML layers

(16)

图 1所示,区域1中主方向是z方向,αx的数值不再是零,而与原衰减函数成比例, 比例系数p(x/z)就是稳定性因子.下面推导多轴近似完全匹配层(M-NPML)

(17)

对于CFS-NPML, 在区域1和2中的处理方法与上面不一样.同样以区域1为例,变换函数如下:

(18)

将式(18)代入式(6)中:

(19)

以上即NPML, CFS-NPML, M-NPML和MCFS-NPML的推导过程和变化方式.在边界区域引入衰减因子来实现对波的吸收,此处选择Groby和Tsogka(2006)改进后的指数型衰减函数,此外引入频移因子、尺度因子和稳定性因子来提升边界的性能.各个因子取值如下:

(20)

L是完全匹配层的厚度, R是理论边界反射系数.l为PML区计算点与内部边界的距离.对于稳定性因子取值应当小于1,大小视介质的复杂程度而改变,不宜太大,否则会导致虚假反射增强.

2 吸收效果及稳定性分析

方程(1)的解有如下形式:

(21)

当进行复拉伸变化之后得到新的波解:

(22)

u0为偏振矢量,kx为波矢量,kx=(cosθ)/cc为相速度,θ为波前面的法向与x轴的夹角.对比式(21)和(22),复拉伸变换之后多出与衰减函数相关项.衰减函数在式(20)中可以看出呈指数增长,波解随着衰减函数积分值的增加而迅速减小,最后趋于0.观察式(22)可知近似完全匹配层存在的两个问题.一是角频率ω在衰减项分母上,当入射波频率很低时产生奇异值,造成不稳定.同时虚假反射的强弱在于积分值的大小,积分值与x值有关,当波的入射角度很大时,近平行于边界层的方向传播,波并不会传入到匹配层的深处,不能被完全吸收.同时未吸收的波在边界中形成耗散波,造成了边界不稳定.当引入频移因子和尺度因子之后,衰减项如下表示:

(23)

在式(23)中,分母上除了角频率还有引入的频移因子,避免了奇异值产生.同时分母上的尺度因子会降低波传播速度(Zhang and Shen, 2010).当波传播速度在边界中因为尺度因子取值而降低的时候,造成波向法线方向弯曲,使波能够传播到边界更深处.

对于稳定性因子,引入的衰减剖面与最初的衰减剖面呈正交分布,能够将传播到边界中的能量进一步吸收,即使是掠入波,也因为平行于边界的衰减剖面而被吸收.完全匹配层中的能量残留造成匹配层不稳定,同时残余在边界中的能量团也会影响边界的吸收效果.

接下来研究NPML的引入对弹性波动方程的影响和稳定性因子在NPML中的作用.为了推导研究整个方程组的稳定性,需要将式(9)和式(17)合并在一起表示.需注意当衰减因子值为0的时候,方程组依旧写为式(1)的形式.现在将其写为矩阵相乘的形式:

(24)

(25)

将式(24)对空间做傅里叶变换并将等式右侧合并为一项:

(26)

(27)

(28)

系数矩阵A与时间无关,式(26)是线性时不变系统.同时矩阵A是非奇异的,因此该方程组有且只有一个平衡点U(t)=0.通过矩阵A的特征值来判断该平衡点的稳定性情况,如果所有的特征值的实部都是负数,那么该平衡点是渐进稳定的.特征值如下表示:

(29)

式(29)中的特征值第二项对应准横波模型,第三项对应准纵波模型.5个特征值实部均为零.以此为基础,研究当我们引入边界之后特征值的变化来了解其稳定性的变化,因此需要将转换方程与波动方程写为式(26)的形式.我们得到如下表示:

(30)

(31)

(32)

(33)

(34)

式(30)是引入近似完全匹配层并做傅里叶变换之后的方程组形式.给出特征矩阵An的同时引入稳定性因子.波场变换再代入到原方程后,特征值变化为式(35):

(35)

当采取近似完全匹配层的波场处理方法处理原方程(未引入衰减函数)时,方程特征值一共变成13个,原方程中的5个不变,额外地增加了8个零特征值.接下来求取矩阵An的特征值:

(36)

式(36)中,系数为

(37)

对比式(35)和(36)可以看出当引入衰减函数之后原式中的9个零特征值减少为5个,其中一个为-αx是负数,其余7个由式(36)中第三项式子求得.其余7个特征值没有办法直接求取,无法了解衰减函数的引入对剩下的特征值的影响.为了预测出这7个特征值在复平面中的变化情况,采取特征值灵敏度分析(Adhikari and Friswell, 2001; Meza-Fajardo and Papageorgiou, 2008).在这里求取特征值对衰减因子的一阶导数.通过观察衰减因子取0的时候的一阶导数分析特征值向复平面正半轴还是负半轴移动.假如当衰减函数从0开始增加的时候,特征值出现在了复平面右侧(特征值的实部大于0),则整个边界是发散不稳定的.

(38)

由式(38)所示,特征值与弹性常数、波传播方向和稳定性因子有关.将式(38)中的波数消去,更方便观察波传播方向对特征值的影响.

(39)

在弹性动力方程中引入微小的衰减函数会使边界不稳定.不稳定性与波的传播方向有关,将本文数值模拟部分的参数代入式(39)中,化简求取得到的特征值灵敏度方程如下:

(40a)

(40b)

(40c)

(40d)

(40e)

(40f)

其中nx=cosθ, θ即是入射角.如公式(40)所示,在弹性各向同性介质中,将特征值灵敏度化为只与传播方向以及稳定性因子有关的方程.没有采取多轴处理的近似完全匹配层稳定性因子为0.此情况下,如式(40a)所示在波沿x轴传播的过程中特征值导数是大于0的,特征值随着小范围衰减函数的引入是增加的,因为是准纵波模型,特征值的实部为0,引入衰减函数之后特征值实部大于0,方程组是不稳定的,一旦引入稳定性因子,特征值导数变为负数,将特征值往复平面左侧偏移,使特征值实部变为负数.式(40b)与(40a)情况类似,式(40c)当稳定性因子取0时,波沿x轴传播的过程中特征值导数取0,稳定性因子引入之后则特征值向左侧偏移,使之对应项由临界稳定变为渐进稳定.式(40d)—(40f)则在波沿z轴传播的过程中产生与此类似情况.从方程组中可以看出衰减函数的引入有可能使特征值朝复平面两侧偏移,有可能如式(36)第二项朝负半轴方向偏移,也可能朝正半轴偏移.然而稳定性因子的引入总是将与之对应的特征值朝负半轴偏移,这也是稳定性因子能够增强其稳定性的原因.如图 2所示,准横波模型、准纵波模型和特征值为0的情况下其一阶导数差异很小,几乎重合,可以清晰看出稳定性因子的引入将整体朝负半轴移动,任意方向的入射波的特征值灵敏度由于稳定性因子的存在变成负数,而且稳定性因子取值越大,一阶导数在负半轴离0越远.具体稳定性因子在实际地震模拟中如何增强其稳定性以及其增强稳定性的效果如何我们通过模拟结果分析.

图 2 特征值灵敏度图 (a)稳定性因子取0;(b)稳定性因子取0.2. Fig. 2 Plots of the eigenvalue sensitivities (a) The stability factor is taken as 0; (b) The stability factor is taken as 0.2.
3 数值模拟分析

第2节从理论上分析叙述了改进后的近似完全匹配层中各个因子在吸收以及增强稳定性方面的作用.接下来从数值模拟结果来进一步解释验证它们的作用以及因子之间的联系.本节采取的模拟参数为:时间步长为1 ms,加载震源是雷克子波,主频是25 Hz.建立的模型大小是5200 m×800 m(包括边界厚度),网格步长x方向和z方向都是10 m.纵波速度3000 m·s-1,横波速度1400 m·s-1,密度为2000 kg·m-3,PML层数L设为10层,R取为10-6.掠入波在入射角接近90°时产生,因此将震源加载于x=2610 m, z=110 m处.在此处产生的波场近似平行于边界传播.

通过模拟快照,反射波振动曲线以及能量衰减曲线分析稳定性因子、频移因子和尺度因子.

图 3来分析稳定性因子对边界的影响.图 3中标示的1、2、3分别对应来自底边界的反射波,来自顶边界的反射波以及来自底边界的反射波在顶边界的二次反射.图 3a第一幅图是NPML下的波场模拟快照,在顶边界中存在着大量的能量团,可以清晰看到这些能量团开始干扰物理区间中的波场.为了观察顶边界中的残留能量情况,调整阈值放大图像,因此也可以看清来自各个边界的反射波.为方便比较,该模型下所有波场快照均采用同一阈值.在NPML中顶边界出现能量团堆积,而底边界没有,因为传播到底边界的入射波入射角较小.在顶边界中由于波沿边界传播,入射角随着偏移距的增加而增大,在离开震源后迅速形成能量“拖尾”现象.近似平行于边界传播的入射波在边界浅层传播,正如第2节公式(22)所示只传播到浅层的入射波无法被完全吸收,其振幅因为衰减因子被衰减,形成耗散波,导致低频奇异值的产生.随着稳定性因子值逐渐增加,边界对能量的吸收效果加强,因为即使沿着边界传播的入射波无法被NPML吸收,但是稳定性因子控制边界从平行方向对其进行吸收,吸收效果视新增衰减剖面的最大因子值而定,因此稳定性因子值越大,边界越“干净”,在边界深处波能量已经完全被吸收.边界中的能量积累使边界稳定性降低,当这些能量被吸收,边界的稳定性也得到提高.但是以顶边界为例,额外引入一个横向的衰减剖面阻挡了波由物理区间进入边界,因此来自于底边界和顶边界的反射波却被大幅度增强.图 3c直观地显示了反射波的强度变化,稳定性因子取0.3时反射波振幅变为原来的5倍.以上的结论在图 3b中得到很好的印证.虚线代表物理区间的能量变化曲线,实线则是整个区间的能量变化曲线,两者之间的差值即代表边界中的能量残留.黑色实箭头代表NPML下边界中残留的能量,而引入稳定性因子之后实线和虚线重合,印证前面所说的边界中能量被进一步吸收.红色实箭头比黑色实箭头的距离大很多,该现象表明M-NPML对边界中能量进一步吸收,防止能量从边界外溢影响物理区间的计算.蓝色箭头所示的物理区间的能量没有进一步减少,表示微弱的反射波没有进入到边界中.

图 3 多轴近似完全匹配层模拟结果 (a)波场快照;(b)能量衰减曲线;(c)反射波振动曲线x=310,z=210. Fig. 3 Simulated results of M-NPML (a) Wave field snapshot; (b) Energy attenuation curve; (c) Reflected wave vibration curve x=310, z=210.

图 4来分析尺度因子与频移因子之间的作用,先看频移因子:在图 4a前三幅图中,随着频移因子的增加,在边界中残留的能量得到衰减,没有出现明显的能量团堆积现象,但是因子值过大在物理区域(图中1和2标示的地方)出现了不稳定现象.图 4d中可以看到频移因子的引入可以削减来自边界的反射波强度,但是在图 4f中却发现在边界中接收到的波的振幅较大.同时观察图 4b,实线与虚线的距离仍然较远,但是在波传播到左右两侧时波能量较NPML有大幅衰减迹象,表明边界处引起的反射能量虚弱了许多,进入边界中的能量增加,但是它们之间的差值也小于NPML,表明入射波在边界中也被进一步吸收.在两侧反射波再次到达两侧边界的时候,吸收能力却被衰弱,同时能量曲线发生震荡现象,边界稳定性减弱,频移因子越大震荡越明显(图中4标识的地方).因此频移因子能够大幅衰减反射波强度,削减边界浅层残留的能量团,但是随着频移因子的增加会产生不稳定现象(与波场快照中1, 2标识处产生的现象对应).

图 4 复频移近似完全匹配层模拟结果 (a)波场快照;(b)(c)能量衰减曲线;(d)(e)反射波振动曲线,x=310,z=210; (f)(g)反射波振动曲线,x=10,z=210. Fig. 4 Simulated results of CFS-NPML (a) Wave field snapshot; (b) (c) Energy attenuation curve; (d) (e) Reflected wave vibration curve, x=310, z=210; (f) (g) Reflected wave vibration curve, x=10, z=210.

分析尺度因子对边界的影响,从波场快照中看出尺度因子有一定削弱大角度入射波产生的能量团的能力,但是尺度因子对边界反射的削弱作用较弱,同时取值过大也会导致不稳定甚至加强反射强度(图 4a中3标识的位置以及图 4e).正如第2节对尺度因子的分析,在图 4g中箭头所示,接收到的波随着尺度因子的增加产生延迟,并且振幅在减小,这也验证之前阐述的尺度因子会降低波在边界中的相速度,值越大速度越低.频移因子增强边界的吸收效果,削弱边界反射,但也会带来一定的不稳定性,尺度因子能够进一步衰减残留在边界深部的能量团,但是过大反而加强了虚假反射.

图 5a左侧两幅图与图 4中标识1所指的地方对比,尺度因子与稳定性因子都减弱了标识1处产生的不稳定现象.图 5a上下图对比,复频移边界下仍然残留的边界深部的能量在MCFS-NPML中被进一步吸收.对于振动曲线,图 5d5e图 5c是一个点接收到的振动曲线,仅仅时间不同,图 5e是将图 5d中矩形中的部分放大.图 5c所示,稳定性因子对反射波的影响在不同的频移因子与尺度因子下是不同的.当尺度因子取值较大时,引入稳定性因子削弱了边界反射同时压制了二次反射波的强烈振动,在频移因子最大值取4的时候反射波却被轻微加强.在进行二次反射波的吸收的过程中,如图 5e所示,MCFS-NPML下的反射波是最弱的,同时在能量衰减曲线中可以看出,当CFS-NPML引入稳定性因子之后,实线与虚线趋于重合,同时由频移因子引起的能量振动现象被完全压制了.

图 5 多轴复频移近似完全匹配层模拟结果 (a)波场快照;(b)能量衰减曲线;(c)(d)(e)反射波振动曲线,x=310,z=210. Fig. 5 Simulated results of MCFS-NPML (a) Wave field snapshot; (b) Energy attenuation curve; (c) (d) (e) Reflected wave vibration curve, x=310, z=210.

下面来观察近似完全匹配层在复杂介质中的模拟情况.

在复杂介质中,边界条件的不稳定现象被放大.在模拟计算的过程中,发现整个系统的不稳定始于边界中,边界中的能量积累并反传回物理区间.在引入频移因子和稳定性因子之后就没有出现该现象,意味着采用近似完全匹配层算法下的系统在没有引入衰减函数的时候是稳定的,即使这时的所有特征值都有0实部,但是引入衰减因子的边界中却不稳定,因为此时特征值出现正数,而引入稳定性因子之后将特征值往负半轴偏移,使系统的稳定性优于引入衰减函数之前.尺度因子对不稳定有一定的压制作用,但能力有限,振动曲线指数增加的转折点仅仅被延后.在能量衰减曲线上可看出频移因子的引入使系统保持稳定,但是在波场快照中可以看到边界中残留的能量仍然十分强烈,稳定性因子引入之后边界中的能量被大幅度吸收,除了边界浅部,能量已经被完全吸收.因此引入稳定性因子对边界的稳定性起着重要作用,但是边界虚假反射却被加强.同时引入尺度因子与频移因子后,削弱了稳定性因子对波的阻挡,加强了对边界深处的残留能量的吸收,进一步增强稳定性的同时也大幅度削弱了虚假反射.

图 6 Marmousi模型下的模拟结果 (a)波场快照(从左到右依次对应图 6b图例中由上到下几种边界条件);(b)能量衰减曲线. Fig. 6 Simulated results under the Marmousi model (a) Wave field snapshots (from left to right corresponding to the top-down boundary conditions in the legend of Fig. 6b); (b) Energy attenuation curve.
4 结论

NPML对波场直接变换,使其具有简单容易实现、不改变控制方程的形式以及对各种介质都有很好的适用性,但因为其对大入射角的地震波吸收不理想导致在正反演的时候依旧采用CPML.为了将NPML更好地应用到正反演中,在该基础上推导出来的MCFS-NPML不仅拥有良好的稳定性,同时增强了对任意方向入射波的吸收能力.经过对方程组波解、特征值灵敏度以及模拟结果进行分析,得到了该匹配层下各个因子之间对完全匹配层的影响,也证明了MCFS-NPML更加稳定以及吸收效果更好.在复杂模拟中只有引入稳定性因子才能清除边界中残留的能量,增强边界的稳定,但是引入稳定性因子的同时吸收效果却下降.为了增强吸收效果引入了尺度因子与频移因子.尺度因子在一定程度上使边界不稳定,频移因子削弱了边界反射,同时避免了低频奇异值的出现,对边界的稳定性有积极贡献,但是当在吸收微弱的二次反射波的时候表现出了不稳定性.尺度因子影响边界中波传播速度,使波向法线方向偏移,使波传播到边界深部.复频移近似完全匹配层对边界浅层的波吸收效果较好,但深部位置的波能量被放大,使边界在某些情况表现出不稳定,因此本文所给出的M-NPML, CFS-NPML在吸收与稳定性上都有自身的缺陷,而MCFS-NPML因为三种因子之间的互补避免了衰减因子引入带来的不稳定以及增强对任意入射方向的波的吸收.

References
Adhikari S, Friswell M I. 2001. Eigenderivative analysis of asymmetric non-conservative systems. International Journal for Numerical Methods in Engineering, 51(6): 709-733. DOI:10.1002/nme.186.abs
Berenger J P. 1994. A perfectly matched layer for the absorption of electromagnetic waves. Journal of Computational Physics, 114(2): 185-200. DOI:10.1006/jcph.1994.1159
Berenger J P. 2004. On the reflection from Cummer's nearly perfectly matched layer. IEEE Microwave and Wireless Components Letters, 14(7): 334-336. DOI:10.1109/LMWC.2004.829272
Chen J Y. 2011. Application of the nearly perfectly matched layer for seismic wave propagation in 2D homogeneous isotropic media. Geophysical Prospecting, 59(4): 662-672. DOI:10.1111/j.1365-2478.2011.00949.x
Chen J Y. 2012. Nearly perfectly matched layer method for seismic wave propagation in poroelastic media. Canadian Journal of Exploration Geophysics, 37(1): 22-27.
Chen K Y. 2010. Study on perfectly matched layer absorbing boundary condition. Geophysical Prospecting for Petroleum (in Chinese), 49(5): 472-477. DOI:10.1360/972010-741
Cheng J W, Mao N B, Lü X C, et al. 2018. Time-domain full waveform inversion with CFS-NPML boundary storage. Oil Geophysical Prospecting (in Chinese), 53(4): 754-764. DOI:10.13810/j.cnki.issn.1000-7210.2018.04.012
Chew W C, Weedon W H. 1994. A 3D perfectly matched medium from modified Maxwell's equations with stretched coordinates. Microwave and Optical Technology Letters, 7(13): 599-604. DOI:10.1002/mop.4650071304
Chew W C, Liu Q H. 1996. Perfectly matched layers for elastodynamics:A new absorbing boundary condition. Journal of Computational Acoustics, 4(4): 341-359. DOI:10.1142/S0218396X96000118
Collino F, Monk P B. 1998. Optimizing the perfectly matched layer. Computer Methods in Applied Mechanics and Engineering, 164(1-2): 157-171. DOI:10.1016/S0045-7825(98)00052-8
Collino F, Tsogka C. 2001. Application of the perfectly matched absorbing layer model to the linear elastodynamic problem in anisotropic heterogeneous media. Geophysics, 66(1): 294-307.
Cummer S A. 2003. A simple, nearly perfectly matched layer for general electromagnetic media. IEEE Microwave and Wireless Components Letters, 13(3): 128-130. DOI:10.1109/LMWC.2003.810124
Feng D S, Yang L Y, Wang X. 2016. The unsplit convolutional perfectly matched layer absorption performance analysis of evanescent wave in GPR FDTD forward modeling. Chinese Journal of Geophysics (in Chinese), 59(12): 4733-4746. DOI:10.6038/cjg20161232
Feng D S, Wang X. 2017. Convolution perfectly matched layer for the finite-element time-domain method modeling of ground penetrating radar. Chinese Journal of Geophysics (in Chinese), 60(1): 413-423. DOI:10.6038/cjg20170134
Festa G, Delavaud E, Vilotte J P. 2005. Interaction between surface waves and absorbing boundaries for wave propagation in geological basins:2D numerical simulations. Geophysical Research Letters, 32(20): L20306. DOI:10.1029/2005GL024091
Groby J P, Tsogka C. 2006. A time domain method for modeling viscoacoustic wave propagation. Journal of Computational Acoustics, 14(2): 201-236. DOI:10.1142/S0218396X06003001
Harari I, Slavutin M, Turkel E. 2000. Analytical and numerical studies of a finite element PML for the Helmholtz equation. Journal of Computational Acoustics, 8(1): 121-137. DOI:10.1142/S0218396X0000008X
Hu F Q. 1996. On absorbing boundary conditions for linearized Euler equations by a perfectly matched layer. Journal of Computational Physics, 129(1): 201-219.
Hu W Y, Abubakar A, Habashy T M. 2007. Application of the nearly perfectly matched layer in acoustic wave modeling. Geophysics, 72(5): SM169-SM175. DOI:10.1190/1.2738553
Huang J P, Huang J Q, Li Z C, et al. 2016. Staggered grid pseudo-spectral numerical simulation of first-order quasi P-wave equation for TTI media. Oil Geophysical Prospecting (in Chinese), 51(3): 487-496. DOI:10.13810/j.cnki.issn.1000-7210.2016.03.010
Komatitsch D, Martin R. 2007. An unsplit convolutional perfectly matched layer improved at grazing incidence for the seismic wave equation. Geophysics, 72(5): SM155-SM167. DOI:10.1190/1.2757586
Kuzuoglu M, Mittra R. 1996. Frequency dependence of the constitutive parameters of causal perfectly matched anisotropic absorbers. IEEE Microwave and Guided Wave Letters, 6(12): 447-449. DOI:10.1109/75.544545
Li Z H, Huang Q H. 2014. Application of the complex frequency shifted perfectly matched layer absorbing boundary conditions in transient electromagnetic method modeling. Chinese Journal of Geophysics (in Chinese), 57(4): 1292-1299. DOI:10.6038/cjg20140426
Lian X M, Shan L Y, Sui Z Q. 2015. An overview of research on perfectly matched layers absorbing boundary condition of seismic forward numerical simulation. Progress in Geophysics (in Chinese), 30(4): 1725-1733.
Liu Y. 2014. Research progress on time-space domain finite difference numerical solution and absorption boundary conditions of wave equation. Oil Geophysical Prospecting (in Chinese), 49(1): 35-46.
Luo Y Q, Liu C. 2018. Absorption effects in nearly perfectly matched layers and damping factor improvement. Oil Geophysical Prospecting (in Chinese), 53(5): 903-913.
Meza-Fajardo K C, Papageorgiou A S. 2008. A nonconvolutional, split-field, perfectly matched layer for wave propagation in isotropic and anisotropic elastic media:stability analysis. Bulletin of the Seismological Society of America, 98(4): 1811-1836. DOI:10.1785/0120070223
Nataf F. 2005. New constructions of perfectly matched layers for the linearized Euler equations. Comptes Rendus Mathematique, 340(10): 775-778. DOI:10.1016/j.crma.2005.04.013
Qin Z, Ren P G, Yao Y, et al. 2009. Improvement of PML absorbing boundary conditions in elastic wave forward modeling. Earth Science-Journal of China University of Geosciences (in Chinese), 34(4): 658-664. DOI:10.3799/dqkx.2009.072
Roden J A, Gedney S D. 2000. Convolution PML (CPML):An efficient FDTD implementation of the CFS-PML for arbitrary media. Microwave and Optical Technology Letters, 27(5): 334-339. DOI:10.1002/1098-2760(20001205)27:5<334::AID-MOP14>3.0.CO;2-A
Teixeira F L, Chew W C. 1998. Extension of the PML absorbing boundary condition to 3D spherical coordinates:scalar case. IEEE Transactions on Magnetics, 34(5): 2680-2683. DOI:10.1109/20.717621
Tian K, Huang J P, Li Z C, et al. 2014. Multi-axial convolution perfectly matched layer absorption boundary condition. Oil Geophysical Prospecting (in Chinese), 49(1): 143-152. DOI:10.1371/journal.pmed.1001323
Wang S D. 2003. Absorbing boundary condition for acoustic wave equation by perfectly matched layer. Oil Geophysical Prospecting (in Chinese), 38(1): 31-34. DOI:10.1007/BF02873153
Xie Z N, Zhang X B. 2017. Weak-form time-domain perfectly matched layer. Chinese Journal of Geophysics (in Chinese), 60(10): 3823-3831. DOI:10.6038/cjg20171012
Xiong Z Q, Mao C Y. 2011. Improved unsplit PML boundary conditions in acoustic wave numerical simulation. Oil Geophysical Prospecting (in Chinese), 46(1): 35-39.
Xue S, Bai D H, Tang J, et al. 2017. Three-dimensional finite difference method for the magnetotelluric secondary field with coupled PML boundary conditions. Chinese Journal of Geophysics (in Chinese), 60(1): 337-348. DOI:10.6038/cjg20170128
Zhang W, Shen Y. 2010. Unsplit complex frequency-shifted PML implementation using auxiliary differential equations for seismic wave modeling. Geophysics, 75(4): T141-T154. DOI:10.1190/1.3463431
Zhang X W, Han L G, Huang L, et al. 2009. A staggered-grid high-order difference method of complex frequency-shifted PML based on recursive integration for elastic wave equation. Chinese Journal of Geophysics (in Chinese), 52(7): 1800-1807. DOI:10.3969/j.issn.0001-5733.2009.07.014
陈可洋. 2010. 完全匹配层吸收边界条件研究. 石油物探, 49(5): 473-477.
成景旺, 毛宁波, 吕晓春, 等. 2018. 非分裂完全匹配层边界存储时间域全波形反演. 石油地球物理勘探, 53(4): 754-764.
冯德山, 杨良勇, 王珣. 2016. 探地雷达FDTD数值模拟中不分裂卷积完全匹配层对倏逝波的吸收效果研究. 地球物理学报, 59(12): 4733-4746. DOI:10.6038/cjg20161232
冯德山, 王珣. 2017. 基于卷积完全匹配层的非规则网格时域有限元探地雷达数值模拟. 地球物理学报, 60(1): 413-423. DOI:10.6038/cjg20170134
黄建平, 黄金强, 李振春, 等. 2016. TTI介质一阶qP波方程交错网格伪谱法正演模拟. 石油地球物理勘探, 51(3): 487-496.
李展辉, 黄清华. 2014. 复频率参数完全匹配层吸收边界在瞬变电磁法正演中的应用. 地球物理学报, 57(4): 1292-1299. DOI:10.6038/cjg20140426
廉西猛, 单联瑜, 隋志强. 2015. 地震正演数值模拟完全匹配层吸收边界条件研究综述. 地球物理学进展, 30(4): 1752-1733.
刘洋. 2014. 波动方程时空域有限差分数值解及吸收边界条件研究进展. 石油地球物理勘探, 49(1): 35-46.
罗玉钦, 刘财. 2018. 近似完全匹配层边界条件吸收效果分析及衰减函数的改进. 石油地球物理勘探, 53(5): 903-913.
秦臻, 任培罡, 姚姚, 等. 2009. 弹性波正演模拟中PML吸收边界条件的改进. 地球科学——中国地质大学学报, 34(4): 658-664.
田坤, 黄建平, 李振春, 等. 2014. 多轴卷积完全匹配层吸收边界条件. 石油地球物理勘探, 49(1): 143-152.
王守东. 2003. 声波方程完全匹配层吸收边界. 石油地球物理勘探, 38(1): 31-34. DOI:10.3321/j.issn:1000-7210.2003.01.007
谢志南, 章旭斌. 2017. 弱形式时域完美匹配层. 地球物理学报, 60(10): 3823-3831. DOI:10.6038/cjg20171012
熊章强, 毛承英. 2011. 声波数值模拟中改进的非分裂式PML边界条件. 石油地球物理勘探, 46(1): 35-39.
薛帅, 白登海, 唐静, 等. 2017. 耦合PML边界条件的三维大地电磁二次场有限差分法. 地球物理学报, 60(1): 337-348. DOI:10.6038/cjg20170128
张显文, 韩立国, 黄玲, 等. 2009. 基于递归积分的复频移PML弹性波方程交错网格高阶差分. 地球物理学报, 52(7): 1800-1807. DOI:10.3969/j.issn.0001-5733.2009.07.014