CN109490947B - 一种高温介质地震波传播模拟方法 - Google Patents
一种高温介质地震波传播模拟方法 Download PDFInfo
- Publication number
- CN109490947B CN109490947B CN201811203247.8A CN201811203247A CN109490947B CN 109490947 B CN109490947 B CN 109490947B CN 201811203247 A CN201811203247 A CN 201811203247A CN 109490947 B CN109490947 B CN 109490947B
- Authority
- CN
- China
- Prior art keywords
- boundary
- displacement
- equation
- physical model
- seismic waves
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 80
- 238000004088 simulation Methods 0.000 title claims abstract description 18
- 238000006073 displacement reaction Methods 0.000 claims abstract description 70
- 239000011435 rock Substances 0.000 claims abstract description 64
- 239000011159 matrix material Substances 0.000 claims description 29
- 239000002245 particle Substances 0.000 claims description 9
- 230000008030 elimination Effects 0.000 claims description 8
- 238000003379 elimination reaction Methods 0.000 claims description 8
- 238000012545 processing Methods 0.000 claims description 8
- 238000009792 diffusion process Methods 0.000 claims description 6
- 238000010521 absorption reaction Methods 0.000 claims description 5
- 239000003637 basic solution Substances 0.000 claims description 5
- 239000000523 sample Substances 0.000 claims description 3
- 230000002159 abnormal effect Effects 0.000 abstract description 3
- 239000000243 solution Substances 0.000 description 13
- 238000004458 analytical method Methods 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 5
- 238000011065 in-situ storage Methods 0.000 description 4
- 230000008859 change Effects 0.000 description 3
- 230000008878 coupling Effects 0.000 description 3
- 238000010168 coupling process Methods 0.000 description 3
- 238000005859 coupling reaction Methods 0.000 description 3
- 239000006185 dispersion Substances 0.000 description 3
- 230000009467 reduction Effects 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- 238000001514 detection method Methods 0.000 description 2
- 230000005489 elastic deformation Effects 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 230000000704 physical effect Effects 0.000 description 2
- 238000006467 substitution reaction Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000005336 cracking Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 239000000463 material Substances 0.000 description 1
- 230000008569 process Effects 0.000 description 1
- 230000002277 temperature effect Effects 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 230000001131 transforming effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/282—Application of seismic models, synthetic seismograms
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明属于地震波模拟领域，具体提供一种高温介质地震波传播模拟方法。本发明旨在解决传统的地震波正演模型为常温介质导致的模拟结果不准确的问题。本发明的模拟方法包括以下步骤：利用边界元法将热弹性波动方程转化为边界元形式方程；利用边界元形式方程求得地震波在均匀岩石物理模型区域的边界的应力值和位移值；将边界的应力值和位移值代入边界元形式方程求得地震波在均匀岩石物理模型内区域的应力值和位移值；根据均匀岩石物理模型区域的边界和内区域的全部应力值和位移值，确定地震波的波形记录。本发明可以准确考察高温对地震波传播的影响，对深层及超深层油气勘探领域拓展提供理论基础，对地球深部观测数据异常现象提供合理的解释。
Description
技术领域
本发明属于地球物理勘探研究领域，具体提供一种高温介质地震波传播模拟方法。
背景技术
地震波波场正演模拟越来越受到地球物理学者们的重视,它的地位也显得越来越重要，正演模拟技术的日趋成熟也使得反演技术在地震勘探中应用的准确性提高。而传统的正演模型往往是常温介质，并没有考虑影响弹性介质物理属性的温度因素，特别是高温状态下介质的弹性属性变化。
众所周知，在不同的温度下，其岩石密度、弹性模量、波速会有较大的差异，尤其对于有裂缝存在的情况下。随着温度的增加，原始的缝洞开裂密度逐渐变大，其必然会引起储层岩性物性参数的变化。所以在实际的深层—超深层地震勘探中，必须考虑到温度的影响，才能更加准确的地反映真实地下构造，故高温介质中地震波传播模拟方法至关重要。
相应地，本领域需要一种高温介质地震波传播模拟方法来解决上述问题。
发明内容
为了解决现有技术中的上述问题，即为了解决传统的地震波正演模型为常温介质导致的模拟结果不准确的问题，本发明提供了一种高温介质地震波传播模拟方法，所述模拟方法包括以下步骤：
步骤S10、利用边界元法将热弹性波动方程转化为边界元形式方程；
步骤S20、利用边界元形式方程求得地震波在均匀岩石物理模型区域的边界的应力值和位移值；
步骤S30、将所述边界的应力值和位移值代入边界元形式方程求得地震波在均匀岩石物理模型内区域的应力值和位移值；
步骤S40、根据均匀岩石物理模型区域的边界和内区域的全部应力值和位移值确定地震波的波形记录。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述步骤S10具体包括以下步骤：
步骤S101、将热弹性波动方程写成频率域形式并转化为二维的情形如下：
公式(1)中，
其中，Δ代表拉普拉斯算子，i为虚数单位，λ和μ代表拉梅常数，u代表质点的位移，ux代表质点在x方向的位移，uz代表质点在z方向的位移，γ代表应力的热力系数，θ代表介质的相对温度，ρ代表介质的密度，ω代表角频率，Xx和Xz分别代表x方向和z方向的体力，i代表第i个离散单元，κ代表热扩散系数，Q代表内热源，k代表导热系数，T0代表介质的初始温度，T代表介质的实时温度，αt代表线膨胀系数，η为中间变量。
步骤S102、将上述热弹性波动方程转化为如下边界元形式方程：
公式(2)中，u和p是公式(2)待求的位移值和应力值，i和j为第i个或者第j个离散单元，N为离散的总单元数，Σ*和V*是热弹性波动方程的应力基本解和位移的基本解，Γ为离散单元的边界。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述步骤S20进一步包括：
步骤S201、给出均匀岩石物理模型参数、观测系统中各检波点及各炮点的位置坐标、边界元法网格单元的各个控制点坐标及区域的划分、震源子波的信息以及均匀岩石物理模型区域边界；
步骤S202、将上述所有信息代入边界元形式方程中，求得地震波在均匀岩石物理模型区域的边界的应力值和位移值。
在上述高温介质地震波传播模拟方法的优选技术方案中，在步骤S202之前，所述模拟方法还包括以下步骤：
利用傅里叶变换将所述的时间域的热弹性波动方程转化为频率域的热弹性波动方程，以便在频率域进行边界元形式方程的计算。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述步骤S202进一步包括：
步骤S2021、利用形函数对所述均匀岩石物理模型的边界进行离散，计算离散点坐标，同时离散外边界吸收层数；
步骤S2022、利用贝塞尔函数计算所述边界元形式方程的积分核函数，得出第一边界矩阵系数；
步骤S2023、利用ASSEM函数对内边界和外边界处理，即将所述第一边界矩阵系数进行组装；
步骤S2024、采用无限元方法进行人工截断边界处理得出第二边界矩阵系数，并将所述第二边界矩阵系数组合到所述第一边界矩阵系数中；
步骤S2025、利用块状高斯消元法对代入步骤S2024中组合后的第一边界矩阵系数的边界元形式方式进行求解，得出地震波在均匀岩石物理模型区域的边界的应力值和位移值。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述步骤S30进一步包括：
步骤S301、将地震波在均匀岩石物理模型区域的边界的应力值和位移值代入边界元形式方程；
步骤S302、利用贝塞尔函数计算上述边界元形式方程的积分核函数，得出第三边界矩阵系数；
步骤S303、利用块状高斯消元法对代入第三边界矩阵系数的边界元形式方式进行求解，得出地震波在均匀岩石物理模型的内区域的应力值和位移值。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述步骤S40进一步包括：利用傅里叶变换将地震波在均匀岩石物理模型的边界的应力值和位移值以及内区域的应力值和位移值从频率域转化为时间域；根据时间域下的地震波在均匀岩石物理模型的边界的应力值和位移值以及内区域的应力值和位移值获得地震波的波形记录。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述均匀岩石物理模型参数包括密度ρ、拉梅常数λ和μ、热扩散系数κ、热膨胀系数αt、单位初始体积常应变比热cε。
在上述高温介质地震波传播模拟方法的优选技术方案中，所述震源子波为高斯子波或者雷克子波，所述震源子波信息包括震源子波的起始频率、中心频率和终止频率。
本领域技术人员能够理解的是，在本发明的技术方案中，高温介质地震波传播模拟方法包括以下步骤：步骤S10、利用边界元法将热弹性波动方程转化为边界元形式方程；步骤S20、利用边界元形式方程求得地震波在均匀岩石物理模型区域的边界的应力值和位移值；步骤S30、将所述边界的应力值和位移值代入边界元形式方程求得地震波在均匀岩石物理模型内区域的应力值和位移值；步骤S40、根据均匀岩石物理模型区域的边界和内区域的全部应力值和位移值确定地震波的波形记录。
通过这样的设置，有如下有益效果：1、由于边界元法的理论基础是边界积分方程理论，其特点是降维、解析与离散相结合、能自动满足无穷远处边界条件，因此该方法只需在区域的边界上划分单元，不必在无穷远边界上划分单元，且计算精度高，适于无限域问题。2、引入温度参量的热弹性波动方程，计及到温度的影响，使结果更加贴切真实情况。本发明实现了基于热弹性波动方程的高温介质地震波传播模拟方法，开展热弹理论研究，可以准确考察高温对地震波传播的影响，对深层及超深层油气勘探领域拓展提供理论基础，对解释地球深部观测数据异常现象提供合理的解释。
附图说明
下面参照附图来描述本发明的优选实施方式，附图中：
图1是本发明一种实施例的高温介质地震波传播模拟方法的主要步骤流程图；
图2是本发明一种实施例的高温介质地震波传播模拟方法的具体步骤流程图；
图3是本发明一种实施例的地震波水平方向的波形图；
图4是本发明一种实施例的地震波垂直方向的波形图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
下面参照附图来描述本发明的优选实施方式。本领域技术人员应当理解的是，这些实施方式仅仅用于解释本发明的技术原理，并非旨在限制本发明的保护范围。
地震波波场正演模拟方法的种类很多，一般将其分为波动方程法、积分方程法、几何射线法三大类。我们基于耦合热弹性理论建立热弹性波动方程，利用边界元数值方法进行高温下地震波传播正演模拟。热弹性耦合理论由热力学基本定律、介质的本构理论和Helmholtz自由能表达式导出的热弹性材料的热传导方程中，除了待定的温度场函数外，还含有应变率。它表明物体上的温度场不仅取决于热源，以及各有关的热力学物性系数和换热边界条件，而且还受到弹性变形应变率的影响，或者说，弹性变形的体积应变率将在一定的程度上改变物体上的热量的传递，因此热传导方程和热弹性运动方程必须耦合求解。
边界元法是对边界积分方程离散求解的现代数值分析方法，属于近似解法之一。边界积分方程从定解问题的控制方程转换而得，其最大的特点就是降低了求解问题的维数，将二维问题化为其边界线上的一维问题，它只以边界变量为基本变量，域内位置可以在需要时根据边界变量求出。基于热弹性理论及边界元法，我们发展一种热介质中弹性波边界元法正演模拟。
如图1所示，图1是本发明一种实施例的高温介质地震波传播模拟方法的主要步骤流程图。参照图1，本发明的高温介质地震波传播模拟方法主要包括以下步骤：
步骤S10、利用边界元法将热弹性波动方程转化为边界元形式方程；
优选地，步骤S10进一步包括：步骤S101、将热弹性波动方程写成频率域形式并转化为二维的情形如下：
公式(1)中，
其中，Δ代表拉普拉斯算子，i为虚数单位，λ和μ代表拉梅常数，u代表质点的位移，ux代表质点在x方向的位移，uz代表质点在z方向的位移，γ代表应力的热力系数，θ代表介质的相对温度，ρ代表介质的密度，ω代表角频率，Xx和Xz分别代表x方向和z方向的体力，i代表第i个离散单元，κ代表热扩散系数，Q代表内热源，k代表导热系数，T0代表介质的初始温度，T代表介质的实时温度，αt代表线膨胀系数，η为中间变量。
如下图所示，在本发明的实施例中常数和一些参数的取值如下：
| ρ(kg·m<sup>-3</sup>) | λ(N·m<sup>-2</sup>) | μ(N·m<sup>-2</sup>) | κ(m<sup>2</sup>/s) | α<sub>t</sub>(1/K) | c<sub>ε</sub>(J/m<sup>3</sup>/K) |
| 2.65 | 4×10<sup>9</sup> | 6×10<sup>9</sup> | 0.089 | 0.33×10<sup>-5</sup> | 1.17×10<sup>2</sup> |
步骤S102、将上述热弹性波动方程转化为如下边界元形式方程：
公式(2)中，Σ*和V*是热弹性波动方程的应力基本解和位移的基本解，u和p是公式(2)待求的位移值和应力值，i和j为第i个或者第j个离散单元，N为离散的总单元数，Γ为离散单元的边界。
一方面，通过引入温度参量的热弹性波动方程，计及到温度的影响，使模拟结果更加贴切实际情况。另一方面，由于边界元法的理论基础是边界积分方程理论，其特点是降维、解析与离散相结合、能自动满足无穷远处边界条件，因此该方法只需在区域的边界上划分单元，不必在无穷远边界上划分单元，且计算精度高，适于无限域问题。
步骤S20、利用边界元形式方程求得地震波在均匀岩石物理模型区域的边界的应力值和位移值。
在本发明的一种实施例中，步骤S20进一步包括：步骤S201、给出均匀岩石物理模型参数、观测系统中各检波点及各炮点的位置坐标、边界元法网格单元的各个控制点坐标及区域的划分、震源子波的信息以及均匀岩石物理模型区域边界。
其中，均匀岩石物理模型参数具体包括：密度ρ、拉梅常数λ和μ、热扩散系数κ、热膨胀系数αt、单位初始体积常应变比热cε共六个参数。这六个参数的信息需要从靶区的原位岩石获取，对靶区的原位岩石进行高温原位热岩石物理实验测量，使求得的参数更加准确、真实可靠。其中，震源子波的信息包括起止频率、中心频率、终止频率等信息。
步骤S202、将上述所有信息代入边界元形式方程中，求得地震波在均匀岩石物理模型区域的边界的应力值和位移值。
需要说明的是，由于在时间域进行边界元形式方程的计算比较复杂，因此为了方便计算，在步骤S202之前，所述方法还包括：利用傅里叶变换将时间域的热弹性波动方程转化为频率域的热弹性波动方程，即公式(1)所示，以便在频率域进行边界元形式方程的计算。
在本发明的一种实施例中，步骤S202进一步包括：
步骤S2021、利用形函数对均匀岩石物理模型的边界进行离散，计算离散点坐标，同时离散外边界吸收层数；
如步骤S201已经给出边界元法网格单元的各个控制点坐标及区域的划分，利用形函数及区域边界的控制点，对区域边界中进行插值，从而计算离散点坐标及离散外边界吸收层数；
步骤S2022、利用贝塞尔函数计算所述边界元形式方程的积分核函数，得出第一边界矩阵系数；
步骤S2023、利用ASSEM函数对内边界和外边界处理，即将所述第一边界矩阵系数进行组装；
步骤S2024、采用无限元方法进行人工截断边界处理得出第二边界矩阵系数，并将所述第二边界矩阵系数组合到所述第一边界矩阵系数中；即在得到第一边界矩阵系数的基础上，再增加扩充矩阵边界系数。
步骤S2025、利用块状高斯消元法对代入步骤S2024中组合后的第一边界矩阵系数的边界元形式方式进行求解，得出地震波在均匀岩石物理模型区域的边界的应力值和位移值。
通过求出地震波在均匀岩石物理模型区域的边界的应力值和位移值，来进行地震波在均匀岩石物理模型内区域的应力值和位移值的求解。
步骤S30、将边界的应力值和位移值代入边界元形式方程求得地震波在均匀岩石物理模型内区域的应力值和位移值；
在本发明的一种实施例中，步骤S30进一步包括：
步骤S301、将地震波在均匀岩石物理模型区域的边界的应力值和位移值代入边界元形式方程；
步骤S302、利用贝塞尔函数计算上述边界元形式方程的积分核函数，得出第三边界矩阵系数；
步骤S303、利用块状高斯消元法对代入第三边界矩阵系数的边界元形式方式进行求解，得出地震波在均匀岩石物理模型的内区域的应力值和位移值。
可以理解的是，要得出多组地震波在均匀岩石物理模型的边界以及内区域的应力值和位移值，需要输入震源子波的某个频率类型的不同频率，得出地震波在震源子波的不同频率下的在均匀岩石物理模型的边界以及内区域的应力值和位移值。如改变给定的震源子波的起始频率：30赫兹、40赫兹、50赫兹等等，计算出地震波在相对应的赫兹下的在均匀岩石物理模型的边界以及内区域的应力值和位移值。
步骤S40、根据均匀岩石物理模型区域的边界和内区域的全部应力值和位移值确定地震波的波形记录。
可以理解的是，由于计算出来的应力值和位移值是在频率下，因此步骤S40进一步包括：利用傅里叶变换将地震波在均匀岩石物理模型的边界的应力值和位移值以及内区域的应力值和位移值从频率域转化为时间域；根据时间域下的地震波在均匀岩石物理模型的边界的应力值和位移值以及内区域的应力值和位移值获得地震波的波形记录。
需要说明的是，上述实施例是在模拟的岩石物理模型为均匀介质下进行的，在岩石物理模型为多层不均匀时，需要引入ANGLE函数对每层地程的倾斜角度进行计算，然后再进行步骤S20、S30和S40。
在一种可能的实施例中，如图2所示，图2是本发明一种实施例的高温介质地震波传播模拟方法的具体步骤流程图。参照图2，高温介质地震波传播模拟方法包括：
1)给定输入的均匀岩石物理模型，其模型包括：密度ρ、拉梅常数λ和μ、热扩散系数κ、热膨胀系数αt、单位初始体积常应变比热cε共六个参数的信息；2)给出观测系统中各检波点及各炮点的位置坐标；3)给出边界元法网格单元的各个控制点坐标及区域的划分；4)给出震源子波的信息，包括起止频率、中心频率、终止频率等信息；5)根据前面4个步骤得到的所有输入信息，输入子模块EDITP中，同时子模块EDITP包含模拟的均匀岩石物理模型区域边界以及模拟震源子波的类型(程序中可选定高斯子波或者雷克子波)；我们选用形函数对均匀岩石物理模型的边界进行离散处理，从而计算边界上的离散点坐标，同时离散外界吸收边界层数；利用傅里叶变换(FFT)对时间域的热弹性波动方程转化为频率域的热弹性波动方程，从而在频率域进行边界元形式方程的离散；调用ASSEM函数对第一边界矩阵系数进行组装，其中第一边界矩阵系数为通过调用贝塞尔函数计算边界元形式方程的积分核函数得出；本程序中采用无限元方法进行人工截断边界处理对第一边界矩阵系数进行重新组装；然后通过块状高斯消元求解，得出地震波在均匀岩石物理模型区域的边界的应力值和位移值；最后将地震波在均匀岩石物理模型区域的边界的应力值和位移值代入边界元形式方程中，通过贝塞尔函数和块状高斯消元求解检波点波场(即地震波在均匀岩石物理模型的区域内的应力值和位移值)；6)将频率域的结果转换为时间域，得出全部节点的位移和应力(即地震波在均匀岩石物理模型的区域内和边界的应力值和位移值)，为后续的反演做准备，从而进行后续的分析。
如图3和图4所示，图3是本发明一种实施例的地震波水平方向的波形图；图4是本发明一种实施例的地震波垂直方向的波形图。参照图3和图4，从特定的输入参数中，我们得到地震波水平和垂直方向的记录，图中我们可以看出，高温介质下的地震波传播记录中一共存在三种波形，即第一纵波、第二纵波(热波)及横波。在不考虑温度效应中，即传统的均匀各项同性介质中的仅仅存在两种记录即横波和纵波记录，我们利用弹性理论模拟地震波在热介质中的传播有新的波形出现，为地震勘探提供新思路。
通过上述描述可以看出，本发明的高温介质地震波传播模拟方法具有以下优点：1、输入的模型岩石物理参数准确，其中密度ρ、拉梅常数λ和μ、热扩散系数κ、热膨胀系数αt、单位初始体积常应变比热cε共六个参数的信息共六个参数的信息需要靶区的原位岩石，从而进行热岩石物理实验测量，求得的参数更加真实可靠。2、利用边界元法，边界元法的理论基础是边界积分方程理论，其特点是降维、解析与离散相结合、能自动满足无穷远处边界条件。该方法只需在区域的边界上划分单元，不必在无穷远边界上划分单元，计算精度高，适于无限域问题。3、引入温度参量的热弹性波动方程，计及到温度的影响，使结果更加贴切真实情况。本发明实现了基于热弹性波动方程的高温介质地震波传播模拟方法，开展热弹理论研究，可以准确考察高温对地震波传播的影响，对深层及超深层油气勘探领域拓展提供理论基础，对解释地球深部观测数据异常现象提供合理的解释。
至此，已经结合附图所示的优选实施方式描述了本发明的技术方案，但是，本领域技术人员容易理解的是，本发明的保护范围显然不局限于这些具体实施方式。在不偏离本发明的原理的前提下，本领域技术人员可以对相关技术特征作出等同的更改或替换，这些更改或替换之后的技术方案都将落入本发明的保护范围之内。
Claims (7)
1.一种高温介质地震波传播模拟方法，其特征在于，包括以下步骤：
步骤S10、利用边界元法将热弹性波动方程转化为边界元形式方程；
步骤S20、利用边界元形式方程求得地震波在均匀岩石物理模型区域的边界的应力值和位移值；
步骤S30、将所述边界的应力值和位移值代入边界元形式方程求得地震波在均匀岩石物理模型内区域的应力值和位移值；
步骤S40、根据均匀岩石物理模型区域的边界和内区域的全部应力值和位移值确定地震波的波形记录；
其中，所述步骤S10具体包括以下步骤：
步骤S101、将热弹性波动方程写成频率域形式并转化为二维的情形如下：
公式(1)中，
其中，Δ代表拉普拉斯算子，i为虚数单位，λ和μ代表拉梅常数，u代表质点的位移，ux代表质点在x方向的位移，uz代表质点在z方向的位移，γ代表应力的热力系数，θ代表介质的相对温度，ρ代表介质的密度，ω代表角频率，Xx和Xz分别代表x方向和z方向的体力，κ代表热扩散系数，Q代表内热源，k代表导热系数，T0代表介质的初始温度，T代表介质的实时温度，αt代表热膨胀系数，η为中间变量，c3为单位初始体积常应变比热；
步骤S102、将上述热弹性波动方程转化为如下边界元形式方程：
公式(2)中，Σ*和V*是热弹性波动方程的应力基本解和位移的基本解，u和p是待求的位移值和应力值，i和j为第i个或者第j个离散单元，N为离散的总单元数，Γ为离散单元的边界；
其中，所述步骤S20进一步包括：
步骤S201、给出均匀岩石物理模型参数、观测系统中各检波点及各炮点的位置坐标、边界元法网格单元的各个控制点坐标及区域的划分、震源子波的信息以及均匀岩石物理模型区域边界；
步骤S202、将步骤S201中给出的所有信息代入边界元形式方程中，求得地震波在均匀岩石物理模型区域的边界的应力值和位移值。
2.根据权利要求1所述的高温介质地震波传播模拟方法，其特征在于，在步骤S202之前，还包括以下步骤：
利用傅里叶变换将时间域的热弹性波动方程转化为频率域的热弹性波动方程，以便在频率域进行边界元形式方程的计算。
3.根据权利要求2所述的高温介质地震波传播模拟方法，其特征在于，所述步骤S202进一步包括：
步骤S2021、利用形函数对所述均匀岩石物理模型的边界进行离散，计算离散点坐标，同时离散外边界吸收层数；
步骤S2022、利用贝塞尔函数计算所述边界元形式方程的积分核函数，得出第一边界矩阵系数；
步骤S2023、利用ASSEM函数对内边界和外边界进行处理，即将所述第一边界矩阵系数进行组装；
步骤S2024、采用无限元方法进行人工截断边界处理得出第二边界矩阵系数，并将所述第二边界矩阵系数组合到所述第一边界矩阵系数中；
步骤S2025、利用块状高斯消元法对代入步骤S2024中组合后的第一边界矩阵系数的边界元形式方程进行求解，得出地震波在均匀岩石物理模型区域的边界的应力值和位移值。
4.根据权利要求3所述的高温介质地震波传播模拟方法，其特征在于，所述步骤S30进一步包括：
步骤S301、将地震波在均匀岩石物理模型区域的边界的应力值和位移值代入边界元形式方程；
步骤S302、利用贝塞尔函数计算上述边界元形式方程的积分核函数，得出第三边界矩阵系数；
步骤S303、利用块状高斯消元法对代入第三边界矩阵系数的边界元形式方程进行求解，得出地震波在均匀岩石物理模型的内区域的应力值和位移值。
5.根据权利要求4所述的高温介质地震波传播模拟方法，其特征在于，所述步骤S40进一步包括：
利用傅里叶变换将地震波在均匀岩石物理模型的边界的应力值和位移值以及内区域的应力值和位移值从频率域转化为时间域；
根据时间域下的地震波在均匀岩石物理模型的边界的应力值和位移值以及内区域的应力值和位移值获得地震波的波形记录。
6.根据权利要求5所述的高温介质地震波传播模拟方法，其特征在于，所述均匀岩石物理模型参数包括介质的密度ρ、拉梅常数λ和μ、热扩散系数κ、热膨胀系数αt、单位初始体积常应变比热cε。
7.根据权利要求6所述的高温介质地震波传播模拟方法，其特征在于，所述震源子波为高斯子波或者雷克子波，
所述震源子波的信息包括震源子波的起始频率、中心频率和终止频率。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811203247.8A CN109490947B (zh) | 2018-10-16 | 2018-10-16 | 一种高温介质地震波传播模拟方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201811203247.8A CN109490947B (zh) | 2018-10-16 | 2018-10-16 | 一种高温介质地震波传播模拟方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109490947A CN109490947A (zh) | 2019-03-19 |
| CN109490947B true CN109490947B (zh) | 2020-01-10 |
Family
ID=65689714
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201811203247.8A Active CN109490947B (zh) | 2018-10-16 | 2018-10-16 | 一种高温介质地震波传播模拟方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109490947B (zh) |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN112329311A (zh) * | 2020-11-10 | 2021-02-05 | 中国石油大学(华东) | 地震波传播有限差分模拟方法、装置及计算机存储介质 |
| CN113468466B (zh) * | 2021-07-23 | 2022-04-15 | 哈尔滨工业大学 | 基于神经网络的一维波动方程求解方法 |
| CN114442154B (zh) * | 2022-04-11 | 2022-06-28 | 中国石油大学(华东) | 变热物性地震波传播模拟方法、系统、设备 |
Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| KR20010074235A (ko) * | 2001-04-12 | 2001-08-04 | 최준성 | 지반-구조물 상호작용의 비선형 시간영역해석을 위한실용적 복합기법 |
| CN104459780A (zh) * | 2014-12-09 | 2015-03-25 | 中国科学院地质与地球物理研究所 | 一种利用波动方程获取地震波路径的方法 |
| CN104597488A (zh) * | 2015-01-21 | 2015-05-06 | 中国石油天然气集团公司 | 非等边长网格波动方程有限差分模板优化设计方法 |
| CN105044771A (zh) * | 2015-08-05 | 2015-11-11 | 北京多分量地震技术研究院 | 基于有限差分法的三维tti双相介质地震波场数值模拟方法 |
| CN108445534A (zh) * | 2018-03-12 | 2018-08-24 | 中海石油(中国)有限公司 | 一种弥散粘滞介质中波场模拟方法 |
-
2018
- 2018-10-16 CN CN201811203247.8A patent/CN109490947B/zh active Active
Patent Citations (5)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| KR20010074235A (ko) * | 2001-04-12 | 2001-08-04 | 최준성 | 지반-구조물 상호작용의 비선형 시간영역해석을 위한실용적 복합기법 |
| CN104459780A (zh) * | 2014-12-09 | 2015-03-25 | 中国科学院地质与地球物理研究所 | 一种利用波动方程获取地震波路径的方法 |
| CN104597488A (zh) * | 2015-01-21 | 2015-05-06 | 中国石油天然气集团公司 | 非等边长网格波动方程有限差分模板优化设计方法 |
| CN105044771A (zh) * | 2015-08-05 | 2015-11-11 | 北京多分量地震技术研究院 | 基于有限差分法的三维tti双相介质地震波场数值模拟方法 |
| CN108445534A (zh) * | 2018-03-12 | 2018-08-24 | 中海石油(中国)有限公司 | 一种弥散粘滞介质中波场模拟方法 |
Non-Patent Citations (3)
| Title |
|---|
| 复杂介质地震波传播模拟中边界元法与有限差分法的比较研究;何彦锋 等;《地球物理学进展》;20130430;第28卷(第2期);第664-678页 * |
| 应用边界元模拟方法分析复杂海底地震散射特征;李绪宣 等;《中国海上油气》;20111231;第23卷(第6期);第357-361、368页 * |
| 弹性波边界元法正演模拟;符力耘 等;《地球物理学报》;19940731;第37卷(第4期);第521-529页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109490947A (zh) | 2019-03-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| You et al. | Spectral element method for dynamic response of transversely isotropic asphalt pavement under impact load | |
| CN109490947B (zh) | 一种高温介质地震波传播模拟方法 | |
| CN106932819B (zh) | 基于各向异性马尔科夫随机域的叠前地震参数反演方法 | |
| Xiang et al. | Experimental investigation of frequency-based multi-damage detection for beams using support vector regression | |
| Ladopoulos | Elastodynamics for Non-linear Seismic Wave Motion in Real-Time Expert Seismology | |
| CN112329311A (zh) | 地震波传播有限差分模拟方法、装置及计算机存储介质 | |
| CN109188519B (zh) | 一种极坐标下的弹性波纵横波速度反演系统及方法 | |
| CN109459787B (zh) | 基于地震槽波全波形反演的煤矿井下构造成像方法及系统 | |
| Zhang et al. | An approach based on level set method for void identification of continuum structure with time-domain dynamic response | |
| CN112949121A (zh) | 一种导波传播特性的求解方法及系统 | |
| Chen et al. | Experiments and computations of solitary wave interaction with fixed, partially submerged, vertical cylinders | |
| Bathe | The finite element method with “overlapping finite elements” | |
| Shen et al. | Dynamic elastic analysis of flexible pavements under moving vehicles: a semi-analytical finite element treatment | |
| Xiao et al. | Dual interpolation boundary face method for 3-D acoustic problems based on binary tree grids | |
| CN111257930B (zh) | 一种黏弹各向异性双相介质区域变网格求解算子 | |
| Zhang et al. | Non-iterative reconstruction of time-domain sound pressure and rapid prediction of large-scale sound field based on IG-DRBEM and POD-RBF | |
| Hernández et al. | Numerical solution of a wave propagation problem along plate structures based on the isogeometric approach | |
| Rajesh et al. | Coupled meshfree and fractal finite element method for mixed mode two‐dimensional crack problems | |
| Zhang et al. | Theoretical study of the dynamic response of a functionally graded hollow cylinder in a half space subjected to plane SH waves | |
| Bause et al. | Transient modeling of ultrasonic guided waves in circular viscoelastic waveguides for inverse material characterization | |
| Sladek et al. | Inverse heat conduction problems in three-dimensional anisotropic functionally graded solids | |
| CN113221409B (zh) | 一种有限元和边界元耦合的声波二维数值模拟方法和装置 | |
| CN106950598B (zh) | 一种偏移速度场可靠性评价方法 | |
| Exadaktylos | A study of the transient fluid flow around a semi-infinite crack | |
| Pepper et al. | A local meshless method for approximating 3D wind fields |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |