张中杰
碳酸盐岩缝洞型储层波场数值模拟
刘 炯 魏修成 陈天胜(中国石化石油勘探开发研究院,北京 100083) 摘 要 缝洞型碳酸盐岩储层作为当今最重要的储层类型之一,一直是油气勘探开发领域的研究重点。但是这类储层往往结构比较复杂,因而地震预测比较困难。本文利用裂缝等效介质理论和随机介质理论建立了在定向裂缝介质中含有随机分布孔洞的各向异性随机模型来描述碳酸盐岩缝洞型储层,并采用伪谱法模拟地震波在模型中的传播。模拟结果显示:由于定向裂缝的作用地震波在缝洞型储层中传播时,地震波会出现横波分裂的现象;此外地震波会因储层中的孔洞发生散射,使得空间波场变得复杂。 关键词 裂缝 孔洞 各向异性 随机介质 Wavefield Simulation in Carbonate Karst Reservoirs LIU Jiong,WEI Xiucheng,CHEN Tiansheng(SINOPEC Exploration & Production Research Institute,Beijing 100083,China) Abstract Carbonate karst reservoir is one of the most important reservoir types in the world,and it has been a research center in oil and gas exploration and development.Because of the complex structure,karst reservoir is difficult to be predicted by traditional seismic exploration technology.Based on fracture equivalent media theory and random media theory,the anisotropic random model is set up to depict carbonate karst reservoir,in which caves and aligned fractures are common.Then pseudospectral method is used to simulate seismic wave propagation in this type of reservoir.Results show that when seismic wave propagate in karst reservoir shear waves will split because of aligned fractures.And on the surfaces of caves,seismic wave will scatter,which makes the wavefield complicated in space. Key words fracture;cave;anisotropy;random media 随着勘探程度的提高,大型构造型油气藏越来越少,勘探目标开始转向复杂储层,如向含裂缝型储层以及裂缝、孔洞结构并存的碳酸盐岩介质转移。 裂缝作为一种复杂的空间结构,大量存在于岩石、地层中。大量的油气勘探实践表明,在储存空间中的裂缝是流体运移的通道,直接关系到油气的产量,同时裂缝在许多储层中也是油气储层的空间,影响储层的油气含量。许多学者对裂缝进行了大量的研究。20世纪80年代,Crampint[1]通过研究发现,地震波在定向裂隙介质中传播时和波在各向异性介质中的传播等效,都会出现快横波和慢横波分裂的现象,并将含定向裂隙的介质称为广泛扩容性各向异性EDA(extensive dilatancy anisotropy)介质。对于一般岩石EDA介质中的众多小裂缝,Hudson[2,3 ]将它们看成是一个个非常扁的椭球体,并用弹性扰动理论推导出裂缝等效各向异性介质的弹性系数与各向同性背景介质的弹性系数、裂缝参数之间的关系,还给出了裂缝中不同充填物对弹性常数的影响。Schoenberg和Sayers[4]将裂缝看成是具有线性滑动边界条件的柔性边界,推导出了裂缝等效各向异性介质的柔性矩阵。随后很多学者运用这两种等效介质理论研究了地震波在裂缝介质中的传播特点。在这些研究中,人们主要从裂缝角度来考虑对地震波的影响,然而实际地层不仅包含裂缝,还可能含有大量尺度不等的其他结构,如孔、洞等。 地层中的孔洞作为区别于背景介质的不均匀结构,往往使得地震波在界面上发生散射、衰减,从而使得波场变得复杂。许多学者对地层中不规则孔洞做了许多工作,研究了它们对地震波传播的影响[5,6]。在以往的孔洞研究中,学者往往用确定性的方法来描述孔洞在空间中的位置,然而对于实际中大量存在的尺度比波长小很多的孔洞,用空间随机分布的方法去表述更为合理。随机模型起源于20世纪60年代。Aki[7]提出了地下介质的随机不均匀性引起的尾波是导致地下振动长久持续的主要原因。Berteusen[8]将随机介质中标量波散射理论应用于大孔径台阵远震P波记录的相位、振幅涨落问题研究。Frankel和Clayton[9]的方法是用数值研究随机弹性介质中纵波的散射衰减。Liu[10]等运用孔隙介质理论和随机理论建立了随机孔隙介质模型,并研究了地震波在其间传播的能量衰减。 前人对单独的裂缝结构和孔洞结构储层已经做了许多研究。然而在实际地层中,地质结构往往不是单一的,如海相碳酸盐岩储层中,由于地质作用,裂缝、孔洞同时存在。为此本文运用Hudson裂缝等效介质理论和随机介质理论建立裂缝介质中含有空间分布小孔洞的各向异性随机介质模型,并采用伪谱法来模拟地震波在该模型中的传播,以此来认识地震波在含微小裂缝和宏观小孔洞的碳酸盐岩储层中的传播规律。 1 各向异性随机介质模型的原理和方法 本部分首先介绍裂缝等效各向异性介质的Hudson理论,然后对弹性波在二维方位各向异性介质中的传播方程做了推导,最后阐述了空间随机介质的建模方法。 1.1 Hudson裂缝等效介质理论 20世纪80年代,Hudson在长波近似、地震波场范围内裂纹位置分布均匀、裂纹在岩石空间中稀疏且彼此不连通的假设前提下,得到了小纵横比扁球体裂缝性质同岩石整体宏观性质之间的关系。在Hudson理论中,含小裂缝的岩石等效弹性常数 可以表示成如下形式: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中: 是各向同性背景介质的弹性常数; 是由于裂缝存在而产生的一阶、二阶修正。 对于垂直裂缝组,裂缝介质显示出横向各向同性的对称性,其总体弹性参数矩阵可以表示为 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 在Hudson理论中,式(2)的各弹性常数的表达式如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中:λ、μ是各向同性背景介质的拉梅常数;ε表示裂缝密度;参数可以表示成 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 垂直裂缝弹性常数表达式中的U1、U3依赖于裂缝内的充填,本文考虑裂缝干燥含气时的状态,表达式如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 1.2 二维方位各向异性介质中地震波的传播 裂缝介质总体弹性系数矩阵(2)是在自身本构坐标系下的表达形式。当裂缝在地质构造作用下发生变化,如层位发生倾斜,即裂缝所在的本构坐标系和观测坐标系存在一定的倾角时,观测坐标系下的弹性矩阵形式会有变化,可以由本构坐标系下的弹性矩阵通过Bond变换获得,其关系如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中:C、C′分别是裂缝介质在本构坐标系和观测坐标系下的弹性系数矩阵;M是坐标转换的Bond矩阵,M′是M的转置。 由式(12)得到裂缝方位各向异性介质的弹性矩阵,其元素一般都不为0。在此基础上,运用弹性介质理论,将应变-位移关系代入关于C′的本构方程中,并将е/еy取为0,可以得到x-z平面内的应力-位移关系如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中:σx、σz分别是介质沿x、z方向的正应力;τyz、τxz、τxy是介质的剪应力;ux、uz是空间质点沿x、z轴的位移分量; 是一般方位各向同性介质的弹性常数。 一般方位各向同性介质中,x-z面内质点的运动方程如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中:ρ表示介质的密度。 方程(13)~(20)构成了方位各向同性介质x-z平面内弹性波传播的控制方程。 1.3 随机介质理论 以往波传播问题的研究多考虑均匀或分层均匀的情况,而实际介质往往是非均匀的。对于实际地下大量存在而且分布不规则的异常介质往往用随机介质模型来描述更接近真实情况。 根据随机过程理论,任意二维空间随机分布量f可以表示成如下平均值和扰动量之和的形式[11]: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中:f0表示f的空间平均值,它是常数;γ(x,z)是在点(x,z)处f相对于平均值的扰动。为了数学上的处理方便,假设空间随机扰动γ(x,z)是均值为0的空间平稳随机过程,即: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 除了均值,人们还往往用方差σ2和自相关函数φ来描述平稳随机过程,它们的表达式如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 上述方程中〈·〉 表示空间平均算子。 根据随机过程理论,γ(x,z)的功率谱就是其自相关函数φ(x,z)的傅立叶变换,所以可以用随机过程的自相关函数用谱展开的方法来构建γ(x,z)空间随机分布。在构建随机介质的过程中自相关函数φ(x,z)的选择有多种,如高斯型、指数型、VonKarman型,本文采用指数型自相关函数,其形式如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 其中a和b分别是随机介质在x和z方向上的自相关长度。 对于均匀背景介质中含孔洞随机分布的模型,可以按以下的方法建立:首先选取指数型函数(25)作为某空间分布量的自相关函数,并选取自相关函数中相关长度,然后用谱展开方法得到函数值的空间随机分布。选取某值作为阈值,当空间某点的分布量大于阈值时认为该点为背景介质区域,否则就是孔洞区域。 2 平面各向异性随机模型的建立和地震波的数值模拟 本节首先建立同时反映裂缝、孔洞性质的各向异性随机模型,并通过地震波波动方程数值模拟的方法,来直观认识地震波在含裂缝、孔洞复杂碳酸盐岩储层中的传播规律。然后建立含碳酸盐岩储层随机各向异性介质的两层地质模型,通过记录地表地震的方法来认识随机各向异性介质模型在地震反射记录上的特点。为了方便认识裂缝、孔洞对地震波的作用,我们对两种模型对应的只含裂缝的情况也做了相应的数值模拟。 2.1 二维各向异性随机模型及其波场模拟 根据前面介绍的Hudson理论建立含微小裂缝的等效各向异性介质,然后以等效的各向异性介质为背景,用谱展开理论建立了含孔洞的各向异性随机介质模型。 为了描述孔洞空间的大小,本文借用孔隙度的概念定义孔洞孔隙度φ如下: 油气成藏理论与勘探开发技术:中国石化石油勘探开发研究院2011年博士后学术论坛文集.4 式中:Vc表示孔洞的体积;Vt是整个模型空间的体积。 2.1.1 二维各向异性随机模型 取自相关函数(25)中的相关尺寸a=b=15m,建立二维各向异性随机模型。如图1所示,背景是含裂缝的等效各向异性介质,其中裂缝主轴水平,偏离观测坐标系x轴-45°。在含裂缝的介质中,背景各向同性介质的密度ρ =2510kg/m3,拉梅常数 λ =7.7351 ×109Pa,μ=23.044×109Pa,裂缝体积密度ε=0.05,其内部假设为干的状态。图1中散布的点代表小尺度孔洞,孔洞中包含物设定为各向同性的泥浆,密度ρc=1580kg/m3,其弹性参数为λc=0.5432×109Pa,μc=1.3665×109Pa,孔洞所占总空间的孔隙度为0.1%。 图1 各向异性随机模型 2.1.2 各向异性随机模型的波场模拟 伪谱法作为波动方程数值模拟的方法之一,主要通过傅立叶变换的方法把物理变量对空间的微分转化为空间频率域中的代数运算,然后再把结果通过傅立叶逆变换转换到物理空间,从而求得对应量的空间微分值[12,13]。理论上其精度可以和有限差分、空间差分精度达到无穷时的情况相当[14]。在本文研究中我们采用伪谱法来模拟地震波在各向异性随机模型中的传播。 首先对图1所示的二维各向异性随机介质空间采用dx=dz=15m的正方形网格进行离散,然后采用伪谱法数值求解。模拟中时间步长取为dt=0.5ms,震源采用心频率取为50Hz的雷克子波,以x方向集中力源的形式安置在模型的中心点。图2是地震波在各向异性随机模型中传播0.325 s时x、y、z3个方向位移分量的波场快照。为了更好地显示孔洞结构对整体波场的影响,我们对只含方位水平裂缝介质的情况也进行了模拟,对应时刻的三分量位移快照如图3所示。 图2 各向异性随机介质中的位移波场快照 从图2、图3的波场快照可以看出准纵波qP、快横波qS1、慢横波qS2的传播,这表明伪谱法可以很好地模拟地震波在各向异性随机介质和裂缝介质中的传播。模拟的结果也表明,由于不同偏振方向的横波传播速度不同,像在裂缝介质中一样,横波在各向异性随机介质中也会产生分裂的现象。与在裂缝介质中不同的是,地震波在各向异性随机模型传播过程中还会出现很多杂乱的散射波。这是因为波遇到孔洞会发生散射,此时每一个小孔洞等效于次生震源,并以此为中心产生与均匀各向异性介质中相似的次生波场。由于模型中含有多个孔洞,多个孔洞产生的次生波场与原波场相互叠加,导致空间总波场变得复杂。 图3 裂缝介质的位移波场快照 2.2 碳酸盐储层各向异性随机模型的地震记录模拟 为了观测缝洞型碳酸盐岩储层在地震记录上的特点,我们设计了含碳酸盐岩储层的两层地质模型,通过伪谱法来模拟地震波在其间的传播,并在地表布置检波器来接收三分量地震位移记录。 2.2.1 含碳酸盐储层的各向异性随机模型 设计含孔洞型碳酸盐岩储层的两层介质模型如图4所示。其中,第一层从0至960m的深度,为碳酸盐岩的各向异性随机介质,其裂缝的对称轴水平,偏离观测坐标系-45°夹角。背景介质的密度ρ1=2400kg/m3,λ1=3.287×109Pa,μ1=17.496×109Pa,裂缝体积密度ε=0.05,孔洞内物质参数同2.2中的各向异性模型,但随机孔洞的孔隙度取为0.05%。第二层从690m至1755m深度,为均匀各向同性介质,其密度ρ2=2650kg/m3,λ2=6.731×109Pa,μ2=32.463×109Pa。 图4 含碳酸盐岩储层的两层地质模型 图5 各向异性随机两层模型的三分量地震记录 2.2.2 碳酸盐储层各向异性随机模型的地震记录模拟 用伪谱法对地震波的传播进行数值模拟。模拟中用dx=dz=15m的正方形网格进行空间离散,时间步长dt=0.5ms;纵波震源采用50Hz的雷克子波,其中心位于x=870m,z=150m的空间点上。在地表每隔15m布置一个检波器接收位移三分量地震记录。整个模拟时间取为0.7s。为了消除人工边界的影响,在模型四周加了完全匹配层吸收边界条件。得到的位移地震记录如图5所示。为了更好地显示各向异性随机模型的反射特征,我们对同样两层模型但第一层是不含孔洞的裂缝各向异性介质的情况也做了模拟,其地震记录如图6所示。 图6 各向异性两层模型的三分量地震记录 从图5、图6的直达准纵波记录qP和直达准横波qS可以看出,虽然在模拟中采用纵波震源激发,但激发的纵波也会产生横波,这是由于在方位裂缝各向异性介质中不同方向的运动是相互耦合的。从模拟结果还可以看出各向异性模型记录中出现了反射准纵波qPqP、反射快横波qPqS1和反射慢横波qPqS2,而碳酸盐岩各向异性随机模型的地震记录上也出现了同样的现象。这是因为准纵波在向下传播时遇到不同介质的分界面,一部分准纵波发生反射并向上传播至地表形成反射准纵波qPqP,一部分准纵波在下行传播过程中遭遇界面时会发生波型转换产生反射横波。由于裂缝结构的作用,反射的横波会分裂成快横波qPqS1和慢横波qPqS2,并先后被地表检波器所接收。从模拟结果还可以看出,在各向异性随机介质的地震记录中出现明显散射现象,而裂缝各向异性介质的地震波记录却没有出现这种现象。这是因为准纵波在传播遭遇到孔洞时,会发生反射和波型转换,产生反射准纵波和反射准横波,这些叠加在裂缝介质的地震记录上,形成各向异性随机介质的复杂地震记录。当孔洞的反射很强时,孔洞的反射波可能会掩盖其他的反射波,如在图5的z方向位移记录上,由于孔洞的影响,反射慢横波变得不明显。 3 结 论 本文运用Hudson裂缝等效介质理论和随机介质理论建立裂缝介质中含有空间分布小孔洞的各向异性随机介质模型来描述海相碳酸盐岩储层中尺度比波长小很多的裂缝和宏观孔洞,采用伪谱法来模拟地震波在该模型中的传播波场;并建立了含碳酸盐岩储层各向异性随机介质的两层模型,模拟了地震波在该模型中的地震反射记录。为了方便认识,同时对不含孔洞的裂缝各向异性模型也做了相应模拟。 两种数值结果都表明,由于介质中定向裂缝的存在,地震波在各向异性随机模型中传播时,会发生横波分裂的现象;地震波会在空间随机分布的孔洞上发生散射,使得整个空间的波场变得复杂,当孔洞的反射很强时,这种散射波可能会在地表记录干扰其他的反射波型,实际解释时要特别加以注意。 参考文献 [1]Crapmpin S.Effective anisotropic elastic constants for wave propagation through cracked solids.Geophysical Journal of the Royal Astronomical,1984,76(1):135~145. [2]Hudson J A.Overall properties of a cracked solid.Mathematical Proceeding of the Cambridge,1980,88:371~384. [3]Hudson J A.Wave speeds and attenuation of elastic waves in material containing cracks.Geophysical Journal of the Royal Astronomical,64(1):133~150. [4]Schoenberg M,Sayers C M.Seismic anisotropy anisotropy of frctred rock.Geophysics,1995,60:204 ~211. [5]Huang B S.A program for two-dimensional seismic wave propagation by the pseudospectrum method.Computers & Geosciences,1992,32:289 ~307. [6]董良国,黄超,刘玉柱,等.溶洞地震反射波特征数值模拟研究.石油物探,2010,49(2):121~124. [7]Aki K.Analysis of the seismic coda of local earthquakes as scattered waves.Journal of Geophysical Research,1969,74:615~631. [8]Berteussen K A,Christoffersson A E S,Husebye E S,and Dahle A.Wave scattering theory in analysis of P wave anomalies at NORSAR and LASA:Geophysical Journal.Royal Astronomical Society,1975,42:403 ~417. [9]Frankel A,Clayton R W.Finite-difference simulations of seismic scattering implications for the propagation of short-period seismic waves in the crust and models of crustal heterogeneity.Journal of Geophysics Research,1986,91,6465 ~6489. [10]Liu J,Ba J,Ma J W,Yang H Z.An analysis of seismic attenuation in random porous media:Science China-Physics,Mechanics& Astronomy,2010,53:628~637. [11]姚姚,奚先.随机介质模型正演模拟及其地震波场分析.石油物探,2002,41(1):31~36. [12]Gazdag J.Modeling of the acoustic wave equation with transform methods.Geophysics,1981,46(6):854~859. [13]Kosloff D D,Baysal E.Forward modeling by a Fourier method.Geophysics,1982,47(10):1402~1412. [14]Forberg B.The pseudospectral method:Comparisons with finite differences for the elastic wave equation.Geophysics,1987,52(4):483~501.
基于碳酸盐岩储层的改进完全匹配层吸收技术
田 坤 李振春 黄建平 李 娜 孔 雪 刘玉金 (中国石油大学(华东)地球科学与技术学院,山东青岛 266555) 基金项目:国家973课题(编号2011CB202402),石油大学创新基金(编号27R1001046A)及(Y090104) 作者简介:田坤,男,在读博士研究生,现从事地震波正演模拟研究。Email:tiankunwudi@163.com。 摘 要:利用波动方程进行数值模拟时,由于介质的计算范围肯定是有限的,这样就会人为地造成 计算边界,不可避免地产生边界反射,因此需要使用吸收边界条件使入射到边界的能量被吸收掉,减小 边界反射所带来的影响。PML吸收层技术已经被证明是非常有效的边界吸收技术,对体波和面波的吸收 都具有非常好的效果,已经被广泛应用于弹性波的数值模拟中。但是在大角度入射即掠射的情况下传统 的PML技术还是存在一定的问题,掠射情况下衰减不够,反射系数比较大,离散后会产生比较严重的 假反射,降低吸收效果,对真实波场产生比较严重的影响,这种影响随入射角、偏移距的增加而增大,而且会对实际应用中的成像、反演等其他处理产生不利影响。而掠射是普遍存在的,比如薄片区域、震 源位置接近研究区域边缘、大偏移距接收等等。本文基于一阶速度-应力方程,提出了一种卷积完全匹 配层(CPML)技术,并采用交错网格有限差分方法对其进行了实现,结果表明CPML技术有效改善了 掠射情况下的吸收效果,尤其是对于散射场等弱能量波场,相对来说效果更是明显。并且在实现过程中 不用分裂变量,应用更加方便简单,易于编程实现,计算卷积时采用递推的形式,不会增加计算量,存 储量也没有太大的变化。 关键词:PML;吸收边界条件;弹性波;数值模拟;卷积完全匹配层技术 An Improved Perfectly Matched Layer Absorbing Boundary Condition for the Carbonate Karst Reservoir Tian Kun,Li Zhenchun,Huang Jianping,Li Na,Kong Xue,Liu Yujin (School of Geoscience,China University of Petroleum(East China),Qingdao 266555,China) Abstract:When we simulate from wave equation,it's can cause computing bounds artificially and the reflections from artificial boundaries inevitably,because the calculation range is limited so it need absorbing boundary conditions to make the energy that incidenting on the bound is absorbed and minimal the impact caused by boundary reflections.The perfectly matched layer(PML)absorbing boundary condition has proven to be a very efficient method for the numerical forward modeling of elastic wave equation to absorb both body waves and surface waves.But at grazing incidence the classical discrete PML method still exist some problems,that is the attenuation is not enough,the reflection coefficient is big to some extent,these canproduce severe spurious reflection after dispersion and reduce the effect to absorb,then cause seriously influences to the real field and these influences are growing with increasing incident angle and offset,besides,these can negatively affect the imaging and inversion and other treatment in practical application.But grazing incidence is natural,for instance in the case of very thin mesh slices,as well as sources located close to the edge of the mesh,or in the case of receivers located at very large offset.In this paper,we present an improved PML absorbing boundary condition at grazing incidence which is called convolution PML(CPML)based on the first-order velocity-stress partial differential wave equation.And it is implemented by a staggered grids finite-difference method.The results show that the CPML technique canobtain a better result in the case of grazing incidence,especially for scattering field and such like this weak reflection wave field.And its implementation is convenient and easy to be used and programmed because of the un-split variable.The cost of the calculation is not increase because of the recursion when convolution is computed.The memory storage is also similar to that of the classical PML Key words:PML;absorbing boundary condition;elastic wave;Forward modeling;CPML 引言 由于计算能力的快速发展,在过去几十年里对地震波在复杂介质中传播的数值模拟方法的研究也越 来越广泛和深入。其中应用最广泛的是有限差分法[1~3],而在有限区域的有限差分正演模拟中,无论是 求解麦克斯韦电磁波方程组,还是求解弹性波动力学方程组,Bérenger(1994)提出的完全匹配层(perfectly matched layer,PML)吸收边界条件是最有效的吸收边界条件[4]。完全匹配层中的波动方程可 以看做是常规的波动方程的推广,波在其中传播时相位改变而振幅随指数衰减;PML和主研究区域弹 性参数相同而衰减系数不同,波阻抗完全匹配,理论上不会发生反射。Chew and Weedon(1994)引入 复数伸展坐标系对PML吸收边界条件进行公式化[5];Rappapport(1995)证明了PML介质等价于在吸 收边界区域引入各向异性介质[6];大量研究[5,7]表明,PML吸收边界条件比指数衰减吸收边界条件[8,9]、 廖氏吸收边界条件[10]、Higdon吸收边界条件[11]和旁轴近似吸收边界条件[12]等具有更优越的吸收性能。Chen et al.(1997)和Wang and Oristaglio(2000)将PML吸收边界条件成功应用到电磁波方程组的求 解中[7,13]。近年来,PML吸收边界条件也被应用到声波和弹性波的有限差分正演模拟中[13~16]。另外,Teixeir and Chew(1999)和Chew and Liu(1996)在柱坐标系和球坐标系中实现了PML吸收边界 条件[17,18]。 但是传统的PML也存在一定的缺陷,离散后的PML层反射系数不严格为0,尤其在大入射角 的情况下更为严重,在这种情况下,相当一部分能量以反射波的形式被送回到主研究区域,而大 入射角的情况是普遍存在的,比如薄片区域、震源位置接近研究区域边缘、大偏移距接收等等。为了克服这个问题,一种修改复坐标变换以改进离散后大入射角吸收效果的PML方法被提了出来 并被应用到麦克斯韦方程组的求解中(Roden and Gedney,2000)并被命名为卷积完全匹配层(CPML)[19]。本文采用这种方法来求解弹性波方程,结果表明CPML技术有效改善了掠射的吸收 效果,并且在实现过程中不用分裂变量,应用更加方便简单,计算卷积时采用递推的形式,不会 增加计算量,存储量也没有太大的变化。 1 CPML基本原理 PML技术本质上是将波动方程在PML层内进行复坐标变换,对于变换后的坐标,方程及其解的形 式是不变的,但是对于原坐标解是衰减的。传统的PML中频率域坐标变换为(以x为例): 国际非常规油气勘探开发(青岛)大会论文集 在CPML中对(3)式进行了扩展,使其形式更一般化: 国际非常规油气勘探开发(青岛)大会论文集 其中αx≥0,κx≥1,可以看出传统PML是CPML在αx=0,κx=1情况下的特例。 将坐标变换关系变换回时域,用 表示1/sx的傅里叶反变换,可以得到 国际非常规油气勘探开发(青岛)大会论文集 根据(4)式可得 国际非常规油气勘探开发(青岛)大会论文集 δ与1和e-atH(t)与1/(a+iω)是傅里叶变换对,所以 国际非常规油气勘探开发(青岛)大会论文集 这样CPML中坐标变换关系就分为两部分,其中第一部分容易计算,只要计算空间偏导除以系数κ 就可以了,最主要的是计算第二部分的卷积。在离散的交错网格下将n时刻的卷积写成 国际非常规油气勘探开发(青岛)大会论文集 交错网格的情况下可以写为 国际非常规油气勘探开发(青岛)大会论文集 其中 国际非常规油气勘探开发(青岛)大会论文集 式中: 国际非常规油气勘探开发(青岛)大会论文集 因为(12)式是简单的指数形式,所以可以把(11)式写成递推的形式: 国际非常规油气勘探开发(青岛)大会论文集 这样在主研究区域进行正常计算,在PML层内通过下面的关系进行坐标变换后的计算就可以了: 国际非常规油气勘探开发(青岛)大会论文集 其中Ψx可以通过(13)式的递推得到。 这就是CPML的基本原理,可以很容易地推广到y,z方向。角上的情形与传统的PML类似,所有 方向同时考虑就可以了。 2 正演模拟计算与对比 为了验证CPML的有效性与优越性,下面对各向同性介质进行了试算,并与传统的PML的计算结 果进行了对比。图1是层状介质模型,采样点为301×301,网格间距5m×5m,上层纵横波速度分别为 3000m/s,1600m/s,下层纵横波速度分别为3300m/s,1900m/s,密度均为2800kg/m3,震源坐标(150,5),道间距5m,采样率0.5ms,采样时间1.5s,图2和图3分别是主频为20Hz的z分量炮记录 和波场快照对比,可以清楚地看到PML在上界面产生了比较强的反射,而且对炮记录有所影响,而 CPML的效果要更好一些。 图4和图5分别是主频为30Hz的z分量炮记录和波场快照对比,也可以看出CPML比PML的 吸收效果要更好一些。图6是3个单道的波形对比图,单道1、2、3坐标分别是(201,5),(251,5),(301,5),图6左侧是0~1.2s的波形对比,右侧是0.3~1.2s的波形对比,也可以 看出CPML要比PML效果好,而且随着偏移距增大,入射角增加,对比越明显,传统的PML中的 反射影响越严重。 图1 层状介质模型 图2 主频为20Hz的炮记录对比 图3 主频为20Hz的0.3s波场快照(上)及局部放大图(下)对比 图4 主频为30Hz的炮记录对比 图7是含散射体介质模型,在(225,15)处有一散射体,散射体纵横波速度分别为2700m/s,1300m/s,其他参数与前相同。图8和图9分别是主频为30Hz的z分量炮记录和波场快照对比,可以看 出PML的虚假反射对散射场有相当的影响,这会对实际应用中的成像、反演等其他处理产生不利影响。 图5 主频为30Hz的0.3s波场快照(上)及局部放大图(下)对比 3 结论 本文基于不分裂的卷积方法在交错网格有限差分情形下对传统的PML技术进行了一般化扩展,推 导了CPML的时域坐标变换关系的递推公式,并通过模型试算进行了对比分析,得到以下几点结论: (1)传统的PML在离散后反射系数不为0,尤其在大入射角、大偏移距的情况下会产生很强的假反射,降低吸收效果,对真实波场产生比较严重地影响,这种影响随入射角、偏移距的增加而增大;(2)改进 后的CPML技术能够增强掠射情况的波场衰减,降低反射系数,有效改善匹配层的吸收效果;(3)对于 含散射体的复杂介质模型,传统的PML也会产生假反射,并更加严重地影响散射场等弱能量波场,而 且会对实际应用中的偏移、反演等处理产生不利影响,而CPML技术同样可以改善吸收效果,降低不利 影响;(4)CPML技术采用不分裂变量的卷积方法进行计算,并通过递推来计算卷积,使其存储量和计算 量与传统PML相比都没有太大变化,不会增加计算成本。 本文的CPML技术能够有效改善掠射的吸收效果,而且不增加计算成本,这对很多如薄片区域、震 源接近研究区域边缘、大偏移距接收等实际情况有很好的应用前景。对于偏移成像等处理中的吸收边界 条件,这也是一种比较好的选择。但是CPML没有改变传统的PML的基本思想,所以它还是存在一些 问题,比如对于一些各向异性介质有固有的不稳定性、离散后匹配层反射系数不严格为零等。本文在匹 配层外围采用的是Dirichlet边界条件,后续研究中可以将其替换为旁轴吸收边界条件以进一步提高离散 后的吸收效果,另外对于二阶位移波动方程的CPML吸收边界条件也有待进一步研究。CPML技术以及 后续相关方法的进一步研究,将有利于西部碳酸盐岩探区复杂近地表速度模型下的正演模拟,尤其是对 碳酸盐岩探区近地表散射波机理认识的研究,可为将来碳酸盐岩探区勘探开发服务。 图6 不同单道的CPML和PML波形对比图 图7 含散射体介质模型 图8 含散射体模型的z分量炮记录对比 图9 含散射体模型的z分量0.3s波场快照(上)及局部放大图(下)对比 参考文献 [1]Alterman Z and Karal F C.Propagation of elastic waves in layered media by finite difference methods.Bulletin of the Seismological Society of America,1968,58:367-398. [2]Madariaga,R..Dynamics of an expanding circular fault.Bulletin of the Seismological Society of America,1976,65: 163-182. [3]Virieux,J.P-SV wave propagation in heterogeneous media:velocity-stress finite-difference method.Geophysics,1986,51: 889-901. [4]Bérenger,J.P .A perfectly matched layer for the absorption of electromagnetic waves.Journal of Computational Physics,1994,114:185-200. [5]Chew W.C.,and W.H.Weedon.A 3-D perfectly matched medium from modified Maxwell's equations with stretched coordinates.Microwave and Optical Technology Letters,1994,7:599-604. [6]Rappaport,C.M..Perfectly matched absorbing boundary conditions based on anisotropic lossy mapping of space.IEEE Microwave Guided Wave Lett.,1995,5:90-92. [7]Chen,Y.H.,Chew,W.C.,Oristaglio,M.L..Application of perfectly matched layers to the transient modeling of subsurface EM problems.Geophysics,1997,62:1730-1736. [8]Marfurt,K.J..Accuracy of finite-difference and finite-element modeling of the scalar and elastic wave equations.Geophysics,1984,49:533-549. [9]Shin,C..Sponge boundary condition for frequency domain modeling.Geophysics,1995,60:1870-1874. [10]Liao,Z.P.,Wong,H.L.,Yang,B.P.,et al..A transmitting boundary for transient wave analysis.Scientia Sinica,1984,27(10):1063-1076. [11]Higdon,R.L.Absorbing boundary condition for elastic waves.Geophysics,1991,56:231-241. [12]Engquist,B.,Majda,A..Absorbing boundary conditions for the numerical simulation of waves.Mathematics of Computation,1977,31(9):629-651. [13]Wang,T.,Oristaglio,M.L.3D simulation of GPR surveys over pipes in dispersive soils.Geophysics,2000,65:1560- 1568. [14]Liu,Q.H.,Tao,J.P..The perfectly matched layer for acoustic waves in absorptive media.Acoust.Soc.Am.,1997,102:2072-2082. [15]Zeng,Y.,He,J.,Liu,Q.H..The application of the perfectly matched layer in numerical modeling of wave propagation in poroelastic media.Geophysics,2001,66:1258-1266. [16]Imhof,M.G..Calculating the seismic effect of 3D underground structures and topography with the finite-difference method.72nd Internat.Mtg.,Soc.Expl.Geophys.,Expanded Abstracts,2002:1939-1942. [17]Teixeira,F.L.,Chew,W.C..On causality and dynamic stability of perfectly matched layers for FDTD simulations.IEEE Transactions on Microwave Theory and Techniques,1999,47:775-785. [18]Chew,W.C.,Liu,Q.H..Perfectly matched layers for elastodynamics:A new absorbing boundary condition.Journal of Computational Acoustics,1996,4:341-359. [19]Roden,J.A.,and S.D.Gedney.Convolution PML(CPML):An efficient FDTD implementation of the CFS-PML for arbitrary media.Microwave and Optical Technology Letters,2000,27:334-339.
科技进步奖前几名有效
科技进步奖在人数之内的都是有效的。 国家科学技术进步奖一等奖单项授奖人数不超过15人,授奖单位不超过10个。 二等奖单项授奖人数不超过10人,授奖单位不超过7个。 特等奖单项授奖人数不超过50人,授奖单位不超过30个。 经评定未授奖的国家自然科学奖、国家技术发明奖和国家科学技术进步奖候选人、候选单位,如果再次以相关项目技术内容推荐须隔一年进行。 授奖条件: 国家科学技术进步奖的奖励范围涉及国民经济的各个行业,是一项覆盖面广泛的科学技术奖。从候选人、候选单位所完成项目的性质来讲,包括了新产品和新技术开发、新技术推广应用、高新技术产业化、企业技术改造及技术进步、技术基础和重大工程建设、重大设备研制中引进消化、吸收国外新技术,或自主开发创新的技术等。 (1)技术创新性突出。 (2)经济效益或者社会效益显著。 (3)推动行业科技进步作用明显。 国家科技进步奖的候选人应当是具备下列条件的项目主要完成人。 (1)提出并确定项目的总体方案。 (2)在解决关键的技术和疑难问题中做出重大技术创新和重要贡献。 (3)在成果转化和推广应用过程做出创造性贡献。 (4)在高技术产业化方面做出重要贡献。候选人按贡献大小排序,并在限额内产生。如果在项目完成中仅从事协调和组织工作的领导,或是从事辅助服务的工作人员,不能作为国家科技进步奖的候选人。 以上内容参考百度百科科学技术进步奖