（固体地球物理学专业论文）双相介质地震波场数值模拟研究.pdf

摘要 双相介质地震波场数值模拟研究 作者简介：李卫志，男，1 9 7 8 年6 月生，2 0 0 4 年9 月从师于李正文教授，于2 0 0 7 年6 月毕业于成都理工大学固体地球物理学专业并获固体地球物理学硕士学位。 摘要 随着地震工程和能源地震勘探的深入发展，为了较精确的描述地下含流体多 孔介质，引入双相介质理论。b i o t ( 1 9 5 6 a 、b ，1 9 6 2 a 、b ) 奠定了双相介质波动传 播理论的基础。b i o t 孔隙介质弹性波理论是建立在宏观水平上的，在具有粘弹性 的流体饱和的孔隙介质中，声波传播的b i o t 线性弹性理论基于以下几点假设：( 1 ) 微细孔隙尺寸远小于地震波波长；( 2 ) 变形很小；( 3 ) 流体相在整个介质中是连续 的，而不连通的孔道可视为固体骨架：( 4 ) 固体骨架是弹性的；( 5 ) 孔隙介质具有 统计上的各向异性，另外还忽略重力的影响。 地震数值模拟是研究地下复杂介质的有效手段，常用的有有限元、有限差分 和虚谱法等，本文重点介绍了有限元和有限差分法的原理、实现方法。在对无边 界介质的地震数值模中，没有引入边界将会影响数值模拟的精度，因此介绍了几 种常用的边界条件。在熟练掌握正演原理及模型建立的基础上，分别建立了若干 符合实际的介质模型进行波动方程的正演模拟。应用实例分析了各正演模拟参数 对模拟的精度及效率的影响，并提出如何在保证正演精度的前提下最大限度的提 高正演效率。 本文通过波场快照分析各复杂介质模型的波动方程响应，进一步加深对复杂 构造波场传播的理解，并在此基础上继续进行双相介质模型的研究。数值模拟结 果表明孔隙流体的存在、固体和流体之间的相互作用会弱化或硬化岩石的力学属 性，从而引起声波或弹性波速度的频散和振幅的衰减，并产生第二类压缩波，裂 缝或裂隙的定向分布、岩层的旋徊性沉积、应力场的定向排列，都会引起传播速 度的各向异性、横波发生分裂等重要现象。 关键词：双相介质，波动方程，有限差分，b i o t 成都理工大学硕士学位论文 s t d u yo nt h en u m e r i c a ls i m u l a t i o no ft h es e i s m i cw a v e f i e l di nt h et w o p h a s e - m e d i u m t h e t w o - p h a s em e d i u mt h e o r y w a si n t r o d u c e di no r d e rt om o r ep r e c i s e d e s c r i p t i o nu n d e r g r o u n dc o n t a i n st h ef l u i dp o r o u sm e d i u mw i t he a r t h q u a k ep r o j e c t a n de n e r g ys e i s m i cs u r v e yt h o r o u g h d e v e l o p m e n t b o i t st h e o r yf o re l a s t i cw a v e si n p o r o u sm e d i aw a se s t a b l i s h e do nm a c r o s c o p i cl e v e l ( b i o t ，1 9 5 6 a ，b ，1 9 6 2 a ，b 1 t h e a n d a s t i ce f f e c t sa r i s ef r o mv i s c o u sj n t e r a c t i o n sb e t w e e nt h ef l u i da n dt h es o i l d m a t r i x t h ef o l l o w i n ga s s u m p t i o n sa r eu s e di nt h et h e o r y ：( 1 ) s e i s m i cw a v e l e n g t hi s l a r g ei nc o m p a r i s o nt ot h ep o r es i z e ( 2 ) t h ed e f o r m a t i o n sa r es m a l l ( 3 ) t h el i q u i d p h a s ei sc o n t i n o u s ，s u c ht h a tp o r ea r ec o n n e c t e da n dt h ed i s c o n n e c t e dp o r ea r ep a r to f t h em a t r i x ( 4 ) t h es o i i dm a t r i xi se l a s t i c ( 5 ) t h em e d i u mi ss t a t i s t i c a l l yi s o t r o p i c ， a n dg r a v i t yf o r c e sa r en e g l e c t e d t h es e i s m i cs i m u l a t i o ni st h er e s e a r c hu n d e r g r o u n dc o m p l e xm e d i u me f f e c t i v e m e t h o d ，t h ef i n i t ee l e m e n t ，f i n i t e d i f f e r e n c e ，a n dp s e u d o s p e c t r a lm e t h o d sa r eu s u a l l y e m p l o y e d i nt h ep a p e r ，t h ef i n i t ee l e m e n ta n dt h ef i n i t e - d i f f e r e n c em e t h o dp r i n c i p l e ， t h er e a l i z a t i o nm e t h o da r ei n t r o d u c e di m p o r t a n t l y i nn u m e r i c a ls i m u l a t i o no fw a v e p r o p a g a t i o nm e d a ，t h ei m p o s i t i o n o fa r t i f i c i a lb o u n d a r i e si n t r o d u c e s p u r i o u s r e f l e c t i o n sw h i c hw i l la f f e c tt h ea c c u r a c yo fn u m e r i a ls u l t i o n s s os e v e r a lb o u n d a r y c o n d i t i o n sa r ei n t r o d u c e d b a s e do nm a s t e r i n gp r i n c i p l e so fw a v ee q u a t i o nf o r w a r d m o d e l i n g ，m a n yc o m p l e xm e d i u mm o d e l ss u i t a b l et ot h er e a l f a c th a v eb e e nb u i l t r e s p e c t i v e l ya n ds u b s e q u e n t l yw a v ee q u a t i o nf o r w a r dm o d e l i n gw a sp r a c t i c e d ，t h e a u t h o ra n a l y z e df o r w a r dm o d e l i n gp a r a m e t e r si n f l u e n c eo nt h em o d e l i n ga c c u r a c y a n de f f i c i e n c yb u ta l s op u tf o r w a r dh o wt oi m p r o v et h ee f f i c i e n c ya tt h em a x i m u m w h i l eg u a r a n t e e i n gf o r w a r dm o d e l i n g a c c u r a c y t h ep a p e ra n a l y z e dw a v ee q u a t i o nr e s p o n s et ot h ec o m p l i c a t e dm e d i u mt h r o u g h i m a g eo fw a v ef i e l d ，t h u st oa c h i e v eab e t t e rc o m p r e h e n do fw a v ep r o p a g a t i o ni n c o m p l i c a t e ds t r u c t u r ea n dc o n t i n u e dt h es u b s e q u e n tp r o c e s s i n ga n dr e s e a r c hf o rt h e t w o - p h a s em e d i u m t h er e s u l tr e v e a l sp o r ef l u i d sa n ds o l i d f l u i di nt e r a c t i o n ( ：l r t s o f t e na n dh a r d e nt h er o c k sm a t r i x t h ev e l o c i t yd i s p e r s i o na n da m p l i t u d ea t t e n u a t i o n o ft h ea c o u s i ca n de l a s t i cw a v e sa n dc a u s et h ec o m p r e s s i o n a lw a v eo ft h es e c o n d k i n d o r i e n t a t i o no ff r a c t u r e so rc r a c k s ，p e r i o d i cs e d i m e n t a r yd e p o s i t so fr o c k ，a n d a l i g n m e n to ft h es t r e s s e sc a u s ev e l o c i t ya n i s o t r o p y ，s h e a r - w a v es p l i t t i n g k e yw o r d s ：t w op h a s e m e d i u m ，w a v ee q u a t i o n ，f i n i t ed i f f e r e n c e ，b i o t 独创性声明 本人声明所呈交的学位论文是本人在导师指导下进行的研究工作及取得的 研究成果。据我所知，除了文中特别加以标注和致谢的地方外，论文中不包含其 他人已经发表或撰写过的研究成果，也不包含为狭得盛叠堡王盔堂或其他教 育机构的学位或证书而使用过的材料。与我一一同工作的同志对本研究所做的任何 贡献均已在论文中作了明确的说明并表示谢意。 签名： 李船 沙u 7 年j 月彳勺e i 学位论文版权使用授权书 本学位论文作者完全了解盛叠堡王太堂有关保留、使用学位论文的规定， 有权保留并向国家有关部门或机构送交论文的复印件和磁盘，允许论文被查阅和 借阅。本人授权盛壑堡王盘堂可以将学位论文的全部或部分内容编入有关数 据库进行检索，可以采用影印、缩印或扫描等复制手段保存、汇编学位论文。 ( 保密的学位沦文在解密后适用本授权书) 签名 瑶崦 ：l o 。7 年，月? o u 第1 章前言 第1 章前言 经典的地震波理论只适合研究固体或流体单相介质中地震波传播规律。然 而，无论砂岩储层还是碳酸盐储层，都是由固体和流体两种部分组成，即由固体 和流体组成的双相介质或多相介质。地震波在双相或多相介质传播时，其传播规 律明显与单相介质中不同。因此，根据双相介质理论，才能可靠准确的确定储层 的厚度、空间展布和储层的孔隙度、渗透率和油气饱和度等储层参数。 1 1 双相介质理论发展现状 1 1 1b i o t 理论发展简介 1 9 5 6 年b i o t 发表了双相介质中地震波传播的基本理论，奠定了双相介质中 弹性波传播规律的基础。b i o t 理论是在g a s s m a n n 方程( 1 9 5 1 ) 出现之后提出的， b i o t 理论要深刻的多，所涉及的面也广泛的多，b i o t 发现了与孔隙介质有关的 慢纵波，以及快纵波、慢纵波及横波的波散和吸收公式。在双相介质中，由于流 体的存在以及固体和流体的相互作用会减弱或硬化岩石的力学属性，弹性波在这 种双相介质中的传播比在单相介质中的传播更为复杂，因此b i o t 双相介质理论 具有经典意义。 b i o t 理论是一种不断发展的理论，从1 9 5 6 年至今，该理论大致可分为三 个重要阶段： 1 9 5 6 年b i o t 发表了流体饱和孔隙固体的弹性波传播理论，首次系统地 提出各向同性饱和多孔隙( 即双相介质) 的弹性波传播理论。1 9 6 2 年b i o t 又 建立了均匀各向异性孔隙线性粘弹性介质中波传播的半唯象理论，也称之为广义 b i o t 理论。b i o t 理论的基本结果是预言了在饱和流体孔隙介质中存在快、慢纵 波和横波，并指出粘滞力控制孔隙流体的相对运动，是弹性波在孔隙介质传播过 程中发生衰减的重要机理【6 1 。 1 9 8 0 年p l o n a 和b e r r y m a n 分别从实验室观测到慢纵波，证实了b i o t 理 论的正确性。但是b i o t 理论的应用范围有一定的局限性，而且不能合理解释高 频散、高衰减的现象。 1 9 9 3 年d v o r k i n 和n u r 针对一维双相系统提出了同时处理b i o t 流动和喷射 流动两种力学机制的统一的b i s q 模型，并将其理论预测结果和实验数据进行比 较，获得了比b i o t 更准确的结果。之后，p a r r a 将这一理论推广到横向各向同性 双相介质情况，杨顶辉等系统地对双相各向异性介质理论中的b i s q 问题进行了 研究。 成都理工大学硕士学位论文 1 1 2 双相介质理论 双相弹性介质理论认为实际的地下介质是由固相、液相组成的，固相的多孔 隙骨架是均匀的、各向同性的弹性固体；液相的充满孔隙空间的物质是具有粘弹 性的、不可压缩的流体。特别是含油储层具有较大的孔隙度，表现出明显的双相 介质性质。 双相介质理论与单相介质理论不同，双相介质理论充分地考虑了介质的结 构、流体与气体的特殊性质、局部特性与整体效应的关系，而纯弹性介质理论不 能详尽地描述岩石对波传播过程的影响，弹性理论只是将岩石简化成单相弹性介 质，将固体和流体、固体和气体、流体和气体，以及它们之间的相互耦合对波的 影响用一些综合岩石参数来描述，这样做在岩石孔隙度很小或孔隙中流体的体压 缩模量和密度很小时是成立的，而当岩石孔隙度较大、孔隙中流体的弹性模量以 及密度均较大时，弹性理论简化就会有偏差。因此，双相介质理论更能准确地描 述实际地层结构和地层性质，自然也就更能适应越来越复杂的油气储藏勘探的实 际需要，从而引起了国内外地震学家和勘探地震学家们的高度重视。双相弹性介 质模型较传统的单相弹性介质假设，更接近实际地层介质，具有可能反映地震信 息与地下岩石的力学性质之间的关系，自然也就更能适应越来越复杂的油气储藏 勘探的实际需要，从而引起了国内外地震学家和勘探地震学家们的高度重视，由 此而发展起来的正演和反演研究具有更好的应用前景。 1 9 5 1 年，g a s s m a n n 提出了关于弹性波在多孔介质中的传播理论，并建立了著 名的g a s s m a n n 方程( 反映了速度与孔隙度之间的定量关系) 。之后，b i o t ( 1 9 5 6 a ， b ) 根据潮湿土壤的电位特性和声学中声波的吸收特性，发展了g a s s m a n n 的流体饱 和多孔隙双相介质理论，奠定了双相介质波动理论的基础【2 - 3 ob i o t 理论充分考 虑了孔隙介质的双相特性，发现了第二类纵波，并指出粘滞力控制孔隙流体的相 对运动，是弹性波在孔隙介质传播过程中发生衰减的重要机理，并于后人的工作 中得到了验证、发展和应用。w h i t e 等人( 1 9 7 5 ) 实验证实了孔隙流体的相对运动是 弹性波衰减的主要机理；p l o n a 实验观测到第二类纵波( 慢纵波) ，证实了b o i t 理论 的正确性。h a m d i 等人用双相介质理论研究了海上沉积物的渗透率；n u t 在综合 前人工作的基础上，较详细地论述了流体饱和度、裂隙密度、孔隙度、孔隙流体 压力与围压、裂隙与孔隙空间之几何形态等等因素对地震波衰减的影响，弹性介 质以及多孔隙岩石内的波传播特点，并指出了双相介质理论在测井储油层评价、 强化回采率、断层检测、圈定地下含蒸汽区域等地震勘探应用前景。n u r 和w a n g 、 w a n g 和n u r ( 1 9 9 2 ) 全面、深入地总结和论述了双相介质中地震波和声波速度的实 验、理论与模型成果；a k b a r 等( 1 9 9 3 ) 研究了p 波衰减与渗透率的变化关系；p a r r a 2 第1 章前言 等( 1 9 9 4 ) 基于b i o t 和均质化理论，研究了波的衰减和频散：s h a r m a ( 1 9 9 6 ) 研究了双 相介质中的面波问题。 在国内，1 1 1 9 0 年代起，有关双相介质中地震波传播理论及其应用的研究也得 到了加强和发展，并已取得了一批重要成果。王尚旭( 1 9 9 0 ) 研究了双相介质中地 震波的传播规律，并实现了双相介质中地震波传播的有限元解法；乔文孝等( 1 9 9 2 ) 研究了声波在两种孔隙介质界面上的反射和透射规律；张应波( 1 9 9 4 ) 探讨了b i o t 理论在地震勘探中的应用；牟永光( 1 9 9 6 ) 通过地震物理模型实验，同时观测并证 实了双相介质中慢速纵波和慢速横波的存在；刘克安等( 1 9 9 7 ) 利用时卷正则化方 法对二维双相介质波动方程中的三个主要参数：孔隙率、固相密度和流相密度同 时进行反演。席道瑛等( 1 9 9 9 ) 研究了在实验室条件下流体的粘性系数与波衰减、 模景和速度色散之间的关系。 1 1 3 双相各向异性研究现状 双相各向异性理论是以双相各向同性理论和地震各向异性理论为基础的，它 产生和发展的重要原因在于这种理论能够更真实地描述地下介质问题”1 。近些 年来，国际上对双相各向异性问题给予了极大的关注，并做了大量的研究工作。 尽管这一问题非常重要，但从双相各向异性理论的发展来看，由于其研究的难度 性和时代性，本世纪8 0 年代末，特别是近几年才是双相各向异性问题研究真正开 始繁荣的时期，而且方兴未艾。 现今对多孔隙双相各向异性理论的研究主要是基于b i o t 模型和等效介质模 型的思想。国际上具有代表性的研究工作有：b i o t ( 1 9 5 7 ，1 9 6 2 a , b ) t 】早期发表的 一系列文章建立的孔隙各向异性介质理论，研究了双相各向异性介质中b i o t 流 动和固体骨架形的力学机制以及地震波传播的规律，这是后来双相各向异性问题 许多研究的基础；1 9 8 8 年k a z i a z o u a 研究了横向各向异性饱和孔隙介质中的格 林函数1 4 4 1 ；p o i m 盱等( 1 9 9 8 ) 利用间接边界元方法模拟了含流体裂隙介质中的散射 波场；s a y e r s ( 1 9 9 9 ) 研究了含流体页岩中依赖于应力的地震各向异性，并指出p 波 的各向异性对渗透率是敏感的】。b e r r y m a n 等( 1 9 8 5 ，1 9 8 6 ，1 9 9 1 ) 提出了几种 等效介质模型i 4 - 删，r a t h o r e 等( 1 9 9 4 ) t s o l t h o m s e n ( 1 9 9 5 ) 【5 ”、h u d s o n 等( 1 9 9 6 ) 5 2 】 综合考虑了弹性固体中的应力场与孔隙、裂缝型介质中的流体场及其相互影响， 进一步发展了孔隙各向异性介质中的等效介质理论。 在国内，牟永光( 1 9 9 6 ) 基于b i o t 理论对双相p t l 、双相e d a 、及双相p t l + e d a 介质中弹性波问题进行了深入研究，给出了双相各向异性介质中弹性波方程的有 限差分方法；刘银斌等( 1 9 9 4 ) 基于广义b i o t 理论，给出了横向各向同性多孑l 介质中 波传播的特征方程；牛滨华等( 1 9 9 4 ) 研究了孔隙各向异性允质中的相速度特性： 成都理工大学硕士学位论文 魏修成( 1 9 9 5 ) 对双相各向异性介质中的波传播特性进行了研究。杨顶辉0 9 9 8 ) 从 微观流场的流体相对流动各向异性速度和b i o t 宏观流动的双相各向异性介质中 弹性波方程出发，实现了双相各向异性介质中关于固体位移和流体“相对流动位 移”的有限元波场模拟。 1 2 数值模拟技术发展现状 数值模拟技术是研究复杂地区地震资料采集、处理和解释的有效辅助手段， 它在石油、天然气、煤、金属等矿产资源以及工程和环境地球物理中得到广泛的 应用。地震数值模拟与地震物理模拟同属于地震正演模拟，即己知地下结构和物 理参数模拟行情况下，预测在地面或井中各观测点所应观测到的地震记录。而地 震勘探的目的则是根据地面或井中各观测点多观测的地震记录来刻画地下介质 结构模型，这是反演过程，却又是建立在正演模型的基础上的。因此地震正演数 值模拟既可以进行正演模拟研究，也是地震反演的基础i l ”。 波动方程数值模拟方法实质上是求解地震波波动方程，因此模拟的地震波场 包含了地震波传播的所有信息。由于波动方程模拟包含了丰富的波动信息，为研 究地震波的传播机理和复杂地层的解释提供了更多的佐证，所以波动方程数值模 拟方法一直在地震模拟中占有重要地位【1 2 】。 全波场波动方程模型正演是种功能很强的波动方程数值模拟方法，在油气 勘探中具有以下的作用【l ”。 1 ) 有助干地球物理方法研究，通过波场快照和正演记录分析来研究各种类 型地震波的传播特征及各种地震响应现象； 2 ) 检验对于油气藏地质构造的认识； 3 ) 在大的勘探项目中核定风险，即通过野外模拟结果来评估和检查旖工的 有效性以减少无效工作； 4 ) 帮助地震资料解释，通过非常清楚和有益的波场模拟来达到去伪存真的 目的； 5 ) 帮助地震资料处理，通过模拟和认识各种波场响应特征来制定合理的处 理流程，选择最佳处理参数。 介质中声波或弹性波场的数值模拟，对于人们理解波动传播规律，解释实际 地震资料，袁征地下介质结构与岩性以及地球资源开发等，均具有重要的理论和 实际意义。随着地震波动理论在天然地震和地震勘探中的应用，地震波理论和计 算机技术的发展，地震数值模拟技术自2 0 世纪6 0 年代以来得到了飞速发展，形 成了目前具有有限差分法、有限元法、虚谱法等各种数值模拟方法的现代地震数 值模拟技术。 4 第1 章前言 有限元法是正演模拟的有效手段，由于剖分的任意性及它所依据的变分原 理，对含有多种介质和自然边界条件的处理，非常方便有效，已成为解决地震波 传播数值模拟的一种重要方法，是目前为止最精确的一种正演模拟方法。l y s m e r 和d r a k e ( 1 9 7 2 ) 最早将有限元法应用到地震数值模拟；s e r o n ( 1 9 9 0 ，1 9 9 6 ) 给 出了弹性波传播有限元模拟方法；p a d o v a n i 等( 1 9 9 4 ) s a r m a 等( 1 9 9 8 ) 给出了 弹性波传播有限元模拟中的无反射条件。 有限差分数值仿真技术是声波或弹性波场数值模拟中最为流行的方法之一， 也是偏微分方程的主要数值解法之一。在各种地震数值模拟方法中，最早出现的 数值模拟方法是有限差分法。a l t e r m a n 和k a m l ( 1 9 6 8 ) 首先将有限差分应用于 层状介质弹性波的数值模拟中；b o o r e ( 1 9 7 2 ) 将有限差分法应用于非均匀介质 地震波传播的模拟；k e l l y 等( 1 9 7 6 ) 研究了用有限差分制作人工合成地震记录 的方法；v i r i e u x ( 1 9 8 6 ) 提出了应用速度一应力一阶方程交错网格有限差分模拟 p s v 波在非均匀介质中的传播。交错网格提高了地震模拟的精度和稳定性，并 消除了部分假象。刘洋等( 1 9 9 8 ) 给出了任意偶数阶精度有限差分法数值模拟， 牟永光、裴正林( 2 0 0 4 ) 对有限差分法等的最新成果进行了总结，编著了三维 复杂介质地震数值模拟等等。有限差分法是一种最常用的正演模拟方法，现已 比较成熟，正向提高精度的方向发展。 伪谱法是一种有效的数值模拟方法，权衡精度和效率，有其他方法不可替代 的优点。在二维介质中，用伪谱法做波传播的正演数值模拟由来已久，由于条件 限制，以前的研究者仅限于二维算法程序的开发和研究。2 0 世纪8 0 年代末， k o s l o f f ( 1 9 8 2 ) 等人用三维声波方程和三维弹性波方程做均匀各向同性介质中波 传播的模拟，并与解析计算结果和超声物理模拟进行比较，证明了方法的正确性。 r e s h e f ( 1 9 8 8 ) 等用伪谱法作了三维声波模拟，所设计的模型在水平和垂直方向 有任意的密度和波速，并应用了吸收边界条件。 1 3 研究思路及研究内容 1 9 5 6 年b i o t 根据潮湿土壤的电位特性和声学中声波的吸收特性，发展 g a s s m a n n 的流体饱和多孔隙介质理论，奠定了双相介质波动理论的基础。b i o t 充分考虑了孔隙介质的双相特性，发现了第二类纵波，指出了粘滞力控制孔隙流 体的相对运动，是弹性波在孔隙介质传播过程中发生衰减的重要机理。之后的一 段时间人们曾对b i o t 理论的真实性有过怀疑，认为不产生第二类纵波是不可能 的，直到p l o n a ( 1 9 8 0 ) 从实验上观测到第二类纵波以后，人们对双相理论的研 究逐渐繁荣起来，直到现在方兴未艾。随着不断的深入研究，对双相介质理论 认识也越来越深刻，也必将对石油地球物理勘探发挥举足轻重的作用。 成都理工大学硬士学位论文 1 3 1 研究思路 在参阅了大量地震波理论进行石油勘探的相关文献后，掌握地震波动方程的 基本理论、地震数值模拟技术和地震波传播规律后，针对不同的地质模型和相应 波动方程，研究对比正演模拟的理论方法，确定数值模拟采用的方法。对数值模 拟方法进行了深入学习和探讨，给出具体的实现步骤，对比不同方法所取得的实 际效果。确定计算理论模型，结合实际测井进行验证，并为反演理论奠定基础。 正演数值模拟技术路线的流程图如下： 1 3 2 研究内容 本文的基本研究内容如下： ( 一) 绪论 介绍了研究双相介质理论和双相各向异性介质的意义、发展现状，对数值模 拟方法的发展现状做了简单的描述。 ( 二) 数值模拟方法的原理 这一部分主要分别描述了限元和有限差分两种数值解法的原理和理论基础。 ( 三) 数值模拟技术的实现 在本章中以弹性波方程为研究对象，分别对有限元法和有限差分的具体实现 步骤，技术路线流程图及程序设计方面做了详细介绍，并几种典型的边界条件做 6 第1 章前言 了一些有意义的讨论。 1 根据均匀介质弹性波方程，采用有限元法以三角单元剖分法，推导了弹 性波有限元方程，并根据有限元的特点介绍了单元矩阵的形成，整体矩阵的组装 方法和存贮方式。形成有限元动力学方程后，根据不同的存贮方式介绍相应的大 型线性方程组的解法，应用有限元单元法进行简单的数值模拟。 2 对均匀介质弹性波方程，应用常规网格差分进行离散，给出对应的差分 格式方程和技术路线流程图。推导均匀介质一阶速度应力方程，采用交错网格 差分进行离散，给出相应的差分格式表达式和技术路线流程图。 3 边界条件在数值模拟是必须要考虑，解决不了边界问题，数值模拟就无 法达到预期的效果，因此边界条件是至关重要的。分别描述吸收边界条件、衰减 边界条件和完全匹配边界条件的方法原理和实现方法，给出每种方法的实现模拟 效果。 ( 四) 复杂地质构造的数值模拟 研究地震波在复杂介质中的传播问题及讨论地震成像方法时，为了方便问题 的求解，经常把时间的固体弹性介质简化为声学介质，只研究纵波的波场特征及 成像方法。这种做法不仅简化了地震波场的构成，同时也对实际问题具有怠好近 似。采用无反射声波方程，采用有限差分数值解法分别对几种模型进行数值模拟。 ( 五) 双相各向同性介质数值模拟 1 9 5 6 年b i o t 发表了流体饱和多孔隙双相介质理论，奠定了双相介质波动理 论的基础。根据b i o t 理论推导双相各向同性纵波方程和一阶速度速度应力方程， 分别利用常规差分和交错网格差分进行数值模拟，分析波场情况。 ( 六) 双相各向异性介质数值模拟 1 9 9 6 年牟永光给出根据b i o t 理论对储层地球物理模型进行了研究，给出了 各向异性介质弹性波方程。根据牟永光给出的有限差分方程进行数值模拟，分析 地震波在双相各向异性介质中的传播情况，建立双相各向异性一阶速度一应力弹 性波方程，利用交错删格差分法进行数值模拟。 1 4 研究成果 本文针对双相介质波动方程的数值解法、地震波场模拟问题以及数值模拟的 边界问题方面展开了一系列的研究。在总结、分析前人研究的基础上，根据导师 的项目课题，对数值模拟方法进行和波动理论进行了深入研究，实现了一套适合 波动方程正演的理论和方法，为观测地震波在地下传播的情况提供了有利的理论 支持和实践指导。主要研究成果体现如下几个方面： 1 本文重点介绍了有限元、有限差分两种数值模拟方法，在有限元法中， 成都理工大学硕士学位论文 整体存贮时，采用了了一种位移一速度解大型动力学方程组的算法，避免了在每 次在时间递推过程中都要解大型线性方程组，大大提高了运行速度。在存贮大型 稀疏矩阵时，为了避免直接存贮占大量内存，使用了一维存贮，最后利用因子化 法解线性方程组。 2 有限差分是数值模拟方法最为常用的方法，也是实现相对简单、模拟效 果较好的一种方法，本文分别对常规有限差分和交错网格差分做了深入研究，重 点介绍了实现的方法和技术路线。在边界的处理上，分别对三种典型的边界条件 进行了描述，并在基础上通过改进，实现任意阶差分的边界处理，模拟效果比较 理想。 本文的创新点是对以往边界处理往往都是针对二阶差分的情况推到了任意 偶数阶，将吸收和衰减边界和混合的方法应用于常规差分的数值模拟中，在交错 网格任意偶数阶差分数值模拟中，分别使用衰减边界和完全匹配边界条件。 3 针对复杂地质构造，采用无反射全声波方程实现任意偶数阶有限差分， 分别对几种典型的模型进行数值模拟，实现自激自收时间剖面的制作，避免了层 问多次波的产生，本方法不但允许速度有纵向变化，而且允许速度有横向变化， 这使得模拟复杂地质构造的地震响应成为可能。 4 根据b i o t 理论，对双相介质理论进行了深入研究，推导了双相各向同性 介质纵波方程和双相各向同性速度一应力弹性波方程，通过数值模拟，研究了地 震波在孔隙双相各向同性介质中的传播情况：在双相介质中除了存在常规的纵 波( 快纵波) 和横波以外，还存在这第二类纵波( 慢纵波) ；慢纵波的速度明 显小于快纵波，而日受耗散系数的影响衰减地很快，所以一般在实际中很难观测 到第二类纵波。 5 双相各向异性介质理论一直是人们关注的热点，因为由此得到的波动方 程更能真实的描述地震波的传播，通过双相各向异性差分方程数值模拟，分析了 地震波的响应情况，通过对耗散系数的调节，观测到慢纵波随着耗散系数的增大 逐渐衰减，通过对弹性参数的调节，得知围体的弹性性质对波的传播有很大的影 响，出现横波分裂等现象。通过交错嘲格差分对双相各向异性介质的波场模拟， 同常规差分的波场模拟进行比较，发现交错网格差分法的精度以及模拟效果要 盘f 。 本章介绍了双相介质发展的现状和数值模拟技术的发展现状，在此基础上提 出了本文的研究内容以及研究思路，最后简要介绍了所取得的研究成果。 第2 章数值模拟方法原理 第2 章数值模拟方法原理 地震数值模拟技术用来模拟波在地层中传播的一种技术，其目的是在给定一 个假设的地层结构以后，预知放置的一系列检波器所记录的地震图。它有助于我 们更精细地描绘地质构造，是地震解释的一种十分有价值的工具，也是地震反演 算法的基础。 随着地震波理论的不断深入研究和计算机技术的飞速发展，数值模拟技术也 得到了快速发展，其中全波动方程方法由于能够完整描述地震波传播规律以及模 拟各种负责介质的地震波传播现象，且对于本构方程、边界条件和震源的形式没 有特定的要求，并允许一定的物性变化而优于积分方程和射线追踪法，成为地震 波数值模拟的主要方法。其中根据数值计算方法的不同，全波动方程方法又可以 为有限元、有限差分、虚谱法等，本章主要介绍有限元和有限差分法的方法原理 和理论基础。 2 1 有限元法 有限元法是结构分析的一种数值计算方法，是矩阵方法在结构力学和弹性力 学等领域中的应用和发展。有限元法借助于矩阵等数学工具，尽管计算工作量很 大，但是整个的分析是一致的，有很强的规律性，因此特别适合于编制计算机程 序来处理。 2 0 世纪5 0 年代中期至2 0 世纪6 0 年代末，有限元法出现并迅猛发展，由于 当时理论尚处于初级阶段，计算机的硬件及软件也无法满足需求，有限元法无法 在工程中得到普遍的应用。从2 0 世纪6 0 年代末2 0 世纪7 0 年代初开始，一些公 司开发出了大型通用的有限元应用程序，这些有限元应用程序以其强大的功能、 简便的操作方法、可靠的计算结果和较高的效率而逐渐形成新的技术商品，成为 结构工程强有力的分析工具。 2 1 1 理论基础 变分法是有限元的数学基础之一，有限元模型最早的数学表达式是建立在变 分法的基础上，变分法在发展单元和解决实际问题方面仍然非常重要，在结构力 学和应力分析领域尤其如此，这两个领域里的近代分析方法几乎是靠有限元法。 变分法的模型通常是要找出一组节点参数值，它使某一特定的称作泛函的积分式 具有驻点值( 即泛函取极大或极小值) 。常常研究的是含多自变量泛函的极小值 问题。 成都理工大学硕士学位论文 m ) = f f “刖昙考蚴一i n n ( 2 1 ) 其中容许函数是q 上具有一阶连续偏导数的函数全体。且在q 的边界r 上满 足某些条件。要使( 2 1 ) 式达到极小的函数，u 必须满足欧拉方程： 一昙( 巧) 一熹( c ) = o ( 2 - 2 ) 其中p ：晏，g ：宴 o x c y 相反的问题式把一个解微分方程的边值问题，化为变分问题。对边值问题可 以设法构造一个对应的泛函，从而化为一个求泛函的极值问题，如下边值问题 等价于变分问题 卜= 蚶留o u 蜊料力m g ”卜曲 ir 0 ：”= 玎 以下列出有限元的计算格式，考虑求解方程 卜；等窘= m 力在区域q ：悱【o ，1 】f 。 ：“| f ：o = 妒( x ) ，“，i , - o = 妒( 工) ，o 蔓xs 1 u ( 0 ，f ) = u ( 1 ，f ) = o ( 2 4 ) ( 2 5 ) 对( 2 4 ) 式使用伽辽金法，假定方程( 2 5 ) 式具有如下形式的近似解 “) = 吼( f ) 苁( z ) ( 2 6 ) ；】 行函数 暑厂 旦砂 芦 旦砂 ， 卜。 叫 n i i a 一缸 弘 声 + (一；誊塑却一 七 “ c = k r 第2 章数值模拟方法原理 苁 ) = 煦，z ， 女一 一1 l ，z 2 峨 ( 2 7 ) 型，x e x k x k 女+ i 一 1 0 ， x 甚 以一1 ，x ；1 】 这里用n + 2 个点把【o ，】分割成n + 1 个子空间阮+ x k ( k = 1 ，2 n + i ) ，这里 所谓的单元或元素，要求这些元素互不重叠。应用( 2 6 ) 式，有 j ( “卜，) 绣级= 0 ，5 = l ，2 ， ( 2 8 ) 即 肿一“一胞p = 0 ( 2 - 9 ) 对变量x 分部积分，并应用边界条件“( 0 ，) = u ( 1 ，) = 0 f 嗔+ “，一胤p = o ( 2 1 0 ) 将( 2 6 ) 式代入( 2 1 0 ) 得 崮( 掣触吼譬豢一能户 协 ls ：l ，2 ，3 n 对每个基本单元 x ix k 】做单元分析，然后在总体合成，最终得到二阶常微分 方程组 m 舀+ k c e = f ( t 1 ( 2 - 1 2 ) 对于( 2 1 2 ) 式可用中心差分求解。 2 1 2 有限元方法的分析过程 1 结构或求解区域的离散化。 有限元方法的思想是化整为零，分散分析，再集零为整。因此，对个结构 进行：有限元分析的第步是将其进行离散，也是根据求解问题的1 i 同精度要求、 效能要求等诸多因素，将整个结构划分有限个单元，单元与边界之间通过结点连 接。因而对于一个具体的有限元分析问题，首先要用适当的有限元把结构进行剖 分，并确定单元的数量、类型、大小和布置。 2 适当的插值模式。 由于在任意给定的约束作用下，问题的准确解为未知，所以我们假设用单元 成都理工大学硕士学位论文 内的一些适当解来近似未知解，从计算的观点看，假设的解必须简单，而且应当 满足一定的收敛性要求。 3 单元分析。 即进行单元刚度矩阵和荷载向量的推导，根据假设的插值模式，利用平衡条 件或适当的变分原理，就可以推导出单元e 的刚度矩阵k 8 和荷载向量f 8 ，形成 单元平衡方程。 k 。p 8 = f 。 4 总体合成。 集合各单元方程可以得到总的平衡方程( 组) 。由于结构是由若干个有限 元组成的，因此，应该把各个单元刚度矩阵和荷载向量按适当的方式进行集合， 从而建立如下形式的总的平衡方程： k p = f 其中，k 为集合刚度矩阵，或称为总体刚度矩阵；p 是整体结构的结点参数 向量，f 是它的节点荷载向量。在不同的领域问题中，p 所代表的物理含义不同， o n ：在固体力学中，p 代表的是结点处的位移，在渗流力学中尸代表的结点处的压 力，热力学中代表的是结点处的温度。 5 引入约束条件。 在总体平衡方程的基础上，按问题的边界条件修改总的平衡方程。考虑了边 界条件，可以把平衡方程表示为 茁p = f 6 方程求解。 对上述方程进行求解，对于线性问可以很容易的解出向量p ；而对非线性问 题则要经过一系列的步骤才能求得解答，每一步都要对矩阵世和荷载向量f 进行 修正。 2 1 3 线性方程组的解法 存建立有限元方程后，最终是落实在对线性方程组的求解。目前解线性方程 组的方法，比较常用的有高斯消元法、三角分解法、迭代法以及共轭梯度法等等。 由于有限元方程在做数值模拟的时候要形成大型的稀疏矩阵，因此在解线性方程 组的时候对其解法有一定的要求，高斯消元一般是解小于1 0 0 0 0 阶的，而晶速度 相对较慢，因此一般不采用来求解有限元方程组。迭代法是一种重要的方法，它 能充分利用系数的稀疏性，减少内存占有量，而且程序实现简单。 利用有限元法求解偏微分方程时，我们一般得到的是如下形式的有限元方程 m p + k p = q ( 2 13 ) 第2 章数值模拟方法原理 对于上面形式的方程，一股情况下是采用e e , c ) 差分逐步积分法，把时l 司区域 【o ，互】分为n 个相等步长址= 五n ，把时刻t 的上式记为 【k 】【尸】f + 【 卅 p 1 = 【q 】i ( 2 1 4 ) 假定在时刻0 ，a t ，2 & t ，t 时刻的解已知，要求t + a t 时刻的解。 把微商用差商代替： 晚2 击 【甩+ 出一垆l + m ( 2 耶 【甩= 面1f 【用l * m - 2 【甩十吼m ( 2 - 1 6 ) 将( 2 1 6 ) 代入( 2 1 4 ) 得 【甩+ 。= 2 【尸】，一【甩一+ 【膨r l q 】，一【网【p 】， a t 2 ( 2 1 7 ) 由此，可以由【p 】。，【p 】，计算出【p 】，。来。 r i o ， 户】o 已给定，【户】o 可由( 2 1 4 ) 式求出，由( 2 - 1 5 ) ，( 2 - 1 6 ) 可得到 【p 】一= 【p 】o 一r 【户】0 + 妻( r ) 2 【j ；】o ( 2 1 8 ) 由( 2 一1 8 ) 可求出【p 】一。，使用【尸1 _ 。，【p 】o 由( 2 1 7 ) 式可以求出【户k 来。 若【q 】_ o ，则有： 旧【尸】+ 妒 = o 可得 e l + & = 2 【尸l 一【p l , 。+ 【 f r 【足】【p l 6 ，2 ( 2 1 9 ) 两边同时乘以【m 】“f ) 2 ，得 器= z 器刊吼一器阢 c 2 圳 为了实现方便 = 器朋= z 等郴， 邝卜器 即可以表示成 【彳】【j p l + 二，= 【口】 p 】， c 】【户l 埘 ( 2 - 2 1 ) 成都理工大学硕士学位论文 f 监1 ：缝! 二幽( 2 - 2 2 ) i o tl址 【只。】= 2 【只卜【一】+ ，2 【 f r l ( 一【k 】【】) 水】+ 愣 础2 m y l ( - 【k h e 1 ) ( 2 - 2 3