CN113435074A - M-ufspml模型、构建方法、智能终端、服务器 - Google Patents

M-ufspml模型、构建方法、智能终端、服务器 Download PDF

Info

Publication number
CN113435074A
CN113435074A CN202110314751.0A CN202110314751A CN113435074A CN 113435074 A CN113435074 A CN 113435074A CN 202110314751 A CN202110314751 A CN 202110314751A CN 113435074 A CN113435074 A CN 113435074A
Authority
CN
China
Prior art keywords
ufspml
order
pml
model
split
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.)
Granted
Application number
CN202110314751.0A
Other languages
English (en)
Other versions
CN113435074B (zh
Inventor
谢志南
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Institute of Engineering Mechanics China Earthquake Administration
Original Assignee
Institute of Engineering Mechanics China Earthquake Administration
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Institute of Engineering Mechanics China Earthquake Administration filed Critical Institute of Engineering Mechanics China Earthquake Administration
Priority to CN202110314751.0A priority Critical patent/CN113435074B/zh
Publication of CN113435074A publication Critical patent/CN113435074A/zh
Application granted granted Critical
Publication of CN113435074B publication Critical patent/CN113435074B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/23Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/11Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Geometry (AREA)
  • Operations Research (AREA)
  • Evolutionary Computation (AREA)
  • Algebra (AREA)
  • Computer Hardware Design (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Geophysics And Detection Of Objects (AREA)

Abstract

本发明属于完美匹配层构建技术领域,公开了一种M‑UFSPML模型、构建方法、智能终端、服务器,对时域二阶线弹性各向异性波动方程的弱形式进行傅里叶变换,得到对应的频域方程;在多轴完全匹配层的概念基础上,采用复坐标延拓技术,对弱形式频域波动方程中的空间变量进行延拓,即得新构建的人工边界,即非分裂卷积多轴完美匹配层。本发明新构建的非分裂卷积多轴完美匹配层比以往的多轴完美匹配层的精度有所提高,同时可保证计算长持时稳定。新构建的非分裂卷积多轴完美匹配层采用非分裂形式,规避了直接离散空间二阶偏导,计算效率更高,其数值结果可保证10万步长以上的数值稳定性,且精度明显优于低阶人工边界。

Description

M-UFSPML模型、构建方法、智能终端、服务器
技术领域
本发明属于完美匹配层构建技术领域,尤其涉及一种M-UFSPML模型、构 建方法、智能终端、服务器,具体涉及一种基于地震各向异性二阶波系统的 M-UFSPML模型构建方法。
背景技术
目前,线性弹性各向异性是地球内部或冰盖(Diez)中非常普遍的属性。它 可能是由各种机制引起的,例如晶体的首选取向,排列的夹杂物或规则的薄层 序列。在各向异性介质中,地震波速度不仅取决于其局部传播和极化方向,还 取决于波的类型(压缩或剪切体波以及Love或Rayleigh面波),还取决于弹性 模型的局部对称性。近年来,数值波模拟作为伴随地震波形层析成像的一部分, 已越来越多地用作研究地球区域尺度各向异性结构的工具。在数值波模拟中, 考虑到各向异性会使数值方案的构建变得复杂,特别是用于无限域截断的技术。 现在众所周知,完全匹配层对于某些种类的各向异性介质会变得不稳定。根据 各向异性弹性参数,从波进入PML区域开始,不稳定性迅速增加。不稳定性与 PML公式的类型无关,例如拆分,未拆分或通过基于速度应力(一阶)或基于 位移(第二阶)的波系统(Komatitsch)的复杂坐标拉伸得出的辅助差分PML 公式。此外,不稳定性也与PML构造中使用的复杂坐标拉伸函数的类型无关。 使用模态分析工具。研究了二维正交各向异性PML波动方程的PML稳定性的 性质,该方程设置为初始值问题。他们发现,当各向异性波方程在计算域中的 慢度曲线不是凸的时,这意味着相速度和群速度的符号相反,所构造的PML将 允许不稳定模式解。而且,他们推导出了一个必要条件和一个单独的充分条件, 均用弹性系数的不等式表示,以排除这种不稳定模式,这限制了PML在各向异 性介质中波传播数值模拟中的应用。一种稳定一阶分裂PML的方法,称为多轴PML(M-PML)。在M-PML中,除了在垂直于PML接口的方向上引入的阻尼 之外,还引入了在与该接口相切的方向上的其他非零阻尼。在他们的研究中, 表明,在2D情况下,较小的切向阻尼水平可使M-PML对于任意各向异性介质 长期稳定。一种数值算法来确定切向阻尼的水平,以稳定2D和3D通用各向异性介质中的一阶分裂PML。尽管通常将切向阻尼设置为比正常阻尼小得多,但 一个重要的缺点是,它使M-PML不再完全匹配。因此,M-PML可解释为标准 PML(无切向阻尼)和具有相同切向和法向阻尼水平的吸收层之间的折衷方案。 为了达到相同的精度水平,M-PML通常在计算效率上更高比标准吸收层基于速 度应力波系统导出的一阶分裂M-PML不能直接用于基于位移的数值模拟方法 中,例如有限元方法,光谱元素法和一些有限差分法。如何利用Newmark时间 步长方案和中点规则之间的等价关系,将一阶分裂PML与基于位移的波动系统 及其有限元/最小二乘离散化在计算域内耦合。根据他们的观察,将他们的一阶 卷积PML扩展到了多轴情况。在他们的UFSPML推导中,使用了复杂的频移 坐标拉伸函数。在他们的一阶M-UFSPML公式中,他们将应力中所有卷积项的 梯度视为附加的体力,这导致相应弱形式方程中存在二阶空间导数。这样的二 阶空间导数不能在基于位移的FE/SE离散化中正确计算,尤其是沿着自由/界面 边界。此外,这种耦合仍然是一个未解决的问题,例如参见最近针对二阶时间 微分方程开发的高阶有效时间步进方案。
为了获得二阶M-PML,直接的方法是将二阶PML扩展到多轴情况。在频 域中开发并实施了二阶UFSPML,并将其扩展到多轴情况。原则上,Li和Bou 的公式不适合分层异质无限域截断,因为同样的原因,他们将应力中所有卷积 项的梯度视为附加的体力。平等。将提出的二阶拆分PML扩展到多轴情况。平 等。分裂位移变量的数量至少是原始变量的五倍。此外,平等。获得了不一致 的结果,与获得的高阶多项式轮廓中具有阻尼参数的M-PML不同,常用的二阶 多项式轮廓中具有阻尼参数的M-PML在长时间仿真中失去了稳定性。
综上所述,现有技术存在的问题是:现有的二阶多轴分离PML(M-SPML) 需要大量的内存和计算量;准确性差。
发明内容
针对现有技术存在的问题,本发明提供了一种M-UFSPML模型、构建方法、 智能终端、服务器。
本发明是这样实现的,一种基于地震各向异性二阶波系统的M-UFSPML模 型构建方法,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法包 括以下步骤:
步骤一,对时域二阶线弹性各向异性波动方程的弱形式进行傅里叶变换, 得到对应的频域方程;
步骤二,在多轴完全匹配层的概念基础上,采用复坐标延拓技术,对弱形 式频域波动方程中的空间变量进行延拓,即得新构建的人工边界。
进一步,所述步骤一中多项式相乘的傅里叶反变换的方法首先通过微调分 母参数以规避极点,并将相乘化为相加的形式,直接通过傅里叶反变换得到其 时域表达式。
进一步,所述步骤二中复坐标延拓函数选用复频移形式。
进一步,所述步骤二中新构建的人工边界被称为非分裂卷积多轴完美匹配 层。
进一步,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法还 包括以下步骤:
第一步,通过对弱形式的基于位移的各向异性波方程进行复杂的拉伸,得 到新型的时域二阶多轴未分裂PML或多轴未分裂频移PML使用经典或频移复 数坐标拉伸函数得出;
第二步,获得的时域二阶多轴未分裂PML或多轴未分裂频移PML以位移 为基础的弱形式,使用有限元或频谱元素方法进行少量修改实现。
进一步,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法的 异性波动方程为:遵循使用的矩阵表示法,各向异性介质由动量方程控制:
Figure BDA0002990702700000041
和应力-应变关系:
σ=C·e (2)
其中u=[ux uy uz]T=[u1 u2 u3]T和f=[f1 f2 f3]T是位移,身体力量向量:
σ=[σxx σyy σzz σyz σxz σxy]T=[σ1 σ2 σy σ4 σ5 σ6]T
e=[σxx eyy ezz eyz exz exy]T=[e1 e2 e3 e4 e5 e6]T
Figure BDA0002990702700000042
在最一般的各向异性情况下,具有三斜对称性的弹性模型,弹性矩阵为:
Figure BDA0002990702700000043
cmn(m,n=1,2,…,6)是由四阶张量cijkl(i,j,k,l=1,2,3)形成的按照以下方式:m对 应于第一对索引i,j,n对应于第二对索引k,l,对应m→i,j和n→k,l如下: 1→1,1;2→2,2;3→3,3;4→2,3;5→1,3;6→1,2;垂直横向各向同性VTI 对称性,水平横向各向同性HTI对称和六角对称各向异性模型的弹性矩阵正交 各向异性对称模型读取:
Figure BDA0002990702700000051
在二维x-y笛卡尔坐标情况下,通过以下方法获得面内波动方程,有 u=[u1 u2]T,σ=[σ1 σ2 σ6]T,e=[e1 e2 e6]T
Figure BDA0002990702700000052
进一步,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法的 通过对PML接受的平面波解的扰动分析,分别获得了必要的和充分的稳定性条 件,为了与(6)中给出的弹性矩阵表示法一致,在x方向上吸收的PML的必 要稳定性条件为:
[(c12+c66)2-c11(c22-c66)]×[(c12+c66)2+c66(c22-c66)]≤0 (7)
足够的稳定性条件表示,如果实现以下两个条件之一,则足以满足PML系 统在x方向上的吸收的稳定性:
(c12+c66)2≤(c11-c66)(c22-c66) (8)
Figure BDA0002990702700000053
本发明的另一目的在于提供一种由所述基于地震各向异性二阶波系统的 M-UFSPML模型构建方法构建的M-UFSPML模型。
本发明的另一目的在于提供一种搭载所述M-UFSPML模型的智能终端。
本发明的另一目的在于提供一种搭载所述M-UFSPML模型的服务器。
综上所述,本发明的优点及积极效果为:本发明提供的基于地震各向异性 二阶波系统的M-UFSPML模型构建方法,构建了适用于各类各向异性介质的基 于二阶位移波动方程的非分裂多轴卷积完美匹配层,该人工边界具有以下优点:
(1)新公式通过复坐标延拓二阶弱形式波动方程,即在延拓过程中统一考 虑了内域和边界条件,保证二者的复坐标延拓的一致性。
(2)新构建的人工边界采用非分裂形式,因此需要更少的计算存储,计算 效率更高。
(3)新构建的人工边界可保证内域和人工边界统一采用有限元或谱元离 散。
(4)新构建的人工边界规避了直接离散空间二阶偏导,可稳定、高精度的 应用至分层各向异性介质中。
针对各向异性非分裂卷积多轴完美匹配层,进行了以下数值实验:
(1)在各向同性半空间均匀弹性介质中,对比了卷积完美匹配层、多轴完 美匹配层以及新构建的非分裂卷积多轴完美匹配层,其结果表明,新构建的非 分裂卷积多轴完美匹配层比以往的多轴完美匹配层的精度有所提高,同时可保 证计算长持时稳定。
(2)分别在全空间均匀分布的第三类、第四类、第五类、锌和磷灰石介质, 以及4层各向异性分层介质,采用新构建的非分裂卷积多轴完美匹配层,其数值 结果可保证10万步长以上的数值稳定性,且精度明显优于低阶人工边界。
本发明通过对弱形式的基于位移的各向异性波系统进行复杂的拉伸,得到 了一种新型的时域二阶多轴未分裂PML(M-UPML)或多轴未分裂频移PML (M-UFSPML)使用经典或频移复数坐标拉伸函数得出。所获得的M-UPML或 M-UFSPML是以位移为基础的弱形式,因此非常适合使用有限元或频谱元素方 法进行少量修改即可实现。新配方具有几个优点。首先,与现有的二阶多轴分 离PML(M-SPML)相比,不再需要二阶和三阶混合微分方程,这大大减少了所 需的内存和计算量。其次,与现有的M-UFSPML相比,不再需要将应力分量中 的卷积项的梯度作为体力以其相应的弱形式公式进行处理,从而避免了潜在的 长期不稳定性和精度损失。使用高阶谱元素离散化和明确的时间步长方案,在 二维波模拟中的理论和实际应用中验证了配方的准确性和通用性。
附图说明
图1是本发明实施例提供的基于地震各向异性二阶波系统的M-UFSPML模 型构建方法流程图。
图2是本发明实施提供的正交各向异性介质III,磷灰石和锌(从上到下) 的慢度曲线(左)和波前(右);灰色曲线与准P波有关,黑色曲线与准S波 有关图。
图3是本发明实施提供的波动时间2s,4s,8s的归一化位移
Figure BDA0002990702700000071
的快照在 各向同性介质中使用不同的无限域截断方法进行的模拟示意图;
图中:(a)9层M-UFSPML;(b)3层UFSPML;(c)9层M-UPML; (d)整个域的能量衰减曲线;(e)在使用不同的无限域截断方法的情况下, 接收端的规范化位移的时程图。
图4是本发明实施提供的对于各向异性介质III中的波模拟示意图;
图中:(a)归一化位移的快照在时间3s
Figure BDA0002990702700000072
(b)归一化位移的快照在 时间6s
Figure BDA0002990702700000073
(c)归一化位移的快照在时间9s的
Figure BDA0002990702700000074
(d)曲线在整个域 中的能量衰减;(e)如果使用了不同的无限域截断方法,则接收机处归一化位 移的两个分量的时程图。
图5是本发明实施提供的对于各向异性介质IV中的波模拟示意图;
图中:(a)归一化位移的快照在时间3s
Figure BDA0002990702700000075
(b)归一化位移的快照在 时间6s
Figure BDA0002990702700000076
(c)归一化位移的快照在时间9s的
Figure BDA0002990702700000077
(d)曲线在整个域 中的能量衰减;(d)曲线在整个域中的能量衰减;(e)如果使用了不同的无 限域截断方法,则接收机处归一化位移的两个分量的时程图。
图6是本发明实施提供的对于各向异性介质V中的波模拟示意图;
图中:(a)归一化位移的快照在时间3s
Figure BDA0002990702700000084
(b)归一化位移的快照在 时间6s
Figure BDA0002990702700000085
(c)归一化位移的快照在时间9s的
Figure BDA0002990702700000086
(d)曲线在整个域 中的能量衰减;(d)曲线在整个域中的能量衰减;(e)如果使用了不同的无 限域截断方法,则接收机处归一化位移的两个分量的时程图。
图7是本发明实施提供的对于各向异性介质锌中的波模拟示意图;
图中:(a)归一化位移的快照在时间60s
Figure BDA0002990702700000081
(b)归一化位移的快照 在时间100s
Figure BDA0002990702700000082
(c)归一化位移的快照在时间200s的
Figure BDA0002990702700000083
(d)整个域 的能量衰减曲线;(e)的时间历程在无限域不同的情况下,接收器的归一化位 移的两个分量使用截断方法图。
图8是本发明实施提供的对于各向异性介质磷灰石中的波模拟示意图;
图8中:(a)归一化位移的快照时间在32s的
Figure BDA0002990702700000087
(b)归一化位移的 快照时间在60s的
Figure BDA0002990702700000088
(c)归一化位移的快照时间在200s的
Figure BDA0002990702700000089
(d) 给定观测点的归一化水平向位移;(e)给定观测点的归一化竖向位移。
图9是本发明实施提供的改进分层模型图。
图10是本发明实施提供的对于分层均质各向异性介质中的波模拟示意图;
图中:(a)归一化位移
Figure BDA00029907027000000810
在时间0.6s;(b)归一化位移
Figure BDA00029907027000000811
在时 间2s;(c)归一化位移
Figure BDA00029907027000000812
在时间4s。
图11是本发明实施提供的整个域中的能量衰减曲线图。
图12是本发明实施提供的接收器上标准化的水平和垂直位移的时间历程如 果使用了不同的无限域截断方法图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白,以下结合实施例, 对本发明进行进一步详细说明。应当理解,此处所描述的具体实施例仅仅用以 解释本发明,并不用于限定本发明。
针对现有技术存在的问题,本发明提供了一种基于地震各向异性二阶波系 统的M-UFSPML模型构建方法,下面结合附图对本发明作详细的描述。
如图1所示,本发明实施例提供的基于地震各向异性二阶波系统的 M-UFSPML模型构建方法包括以下步骤:
S101,对时域二阶线弹性各向异性波动方程的弱形式进行傅里叶变换,得 到对应的频域方程。
S102,在多轴完全匹配层的概念基础上,采用复坐标延拓技术,对弱形式 频域波动方程中的空间变量进行延拓,即得新构建的人工边界。
本发明实施例提供的S101中,所述多项式相乘的傅里叶反变换的方法如下: 首先通过微调分母参数以规避极点,并将相乘化为相加的形式,从而可直接通 过傅里叶反变换得到其时域表达式。
本发明实施例提供的S102中,所述复坐标延拓函数选用复频移形式。
本发明实施例提供的S102中,所述新构建的人工边界被称为非分裂卷积多 轴完美匹配层。
下面结合附图对本发明的技术方案作进一步的描述。
本发明实施例提供的基于地震各向异性二阶波系统的M-UFSPML模型构建 方法通过对弱形式的基于位移的各向异性波方程进行复杂的拉伸,得到了一种 新型的时域二阶多轴未分裂PML(M-UPML)或多轴未分裂频移PML (M-UFSPML)使用经典或频移复数坐标拉伸函数得出;获得的M-UPML或 M-UFSPML是以位移为基础的弱形式,使用有限元或频谱元素方法进行少量修 改即可实现。
1、二阶波方程及PML在各向异性介质中的应用局限性
1.1二阶弹性各向异性波动方程
遵循使用的矩阵表示法,各向异性介质由动量方程控制
Figure BDA0002990702700000101
和应力-应变关系:
σ=C·e, (2)
其中u=[ux uy uz]T=[u1 u2 u3]T和f=[f1 f2 f3]T是位移,身体力量向 量:
σ=[σxx σyy σzz σyz σxz σxy]T=[σ1 σ2 σ3 σ4 σ5 σ6]T
e=[exx eyy ezz eyz exz exy]T=[e1 e2 e3 e4 e5 e6]T
Figure BDA0002990702700000102
在最一般的各向异性情况下,即具有三斜对称性的弹性模型,弹性矩阵为:
Figure BDA0002990702700000103
cmn(m,n=1,2,…,6)是由经典的四阶张量cijkl(i,j,k,l=1,2,3)形成的按照以下方式:m对应于第一对索引i,j,n对应于第二对索引k,l,对应m→i,j和n→k,l 如下:1→1,1;2→2,2;3→3,3;4→2,3;5→1,3;6→1,2。让本发明将自身 限制为正交的各向异性模型对称性,在最近的地震学和地震勘探中起着最重要 的作用。它涵盖了VTI(垂直横向各向同性)对称性,HTI(水平横向各向同性) 对称和六角对称各向异性模型。的弹性矩阵正交各向异性对称模型读取:
Figure BDA0002990702700000111
在二维情况下(x-y笛卡尔坐标),可以通过以下方法获得面内波动方程: 只需省略与z坐标相关的变量和空间变化。然后,本发明有 u=[u1 u2]T,σ=[σ1 σ2 σ6]T,e=[e1 e2e6]T
Figure BDA0002990702700000112
1.2 PML在各向异性介质中应用的局限性
二维平面内一阶分裂PML的稳定性各向异性波模拟。将具有恒定阻尼参数 的PML中的波动视为初始值问题。通过对PML接受的平面波解的扰动分析, 他们分别获得了必要的和充分的稳定性条件,该条件用各向异性模型的弹性参 数的不等式表示。为了与(6)中给出的弹性矩阵表示法一致,在x方向上吸收 的PML的必要稳定性条件为:
[(c12+c66)2-c11(c22-c66)]×[(c12+c66)2+c66(c22-c66)]≤0. (7)
足够的稳定性条件表示,如果实现以下两个条件之一,则足以满足PML系 统在x方向上的吸收的稳定性:
(c12+c66)2≤(c11-c66)(c22-c66). (8)
Figure BDA0002990702700000113
尽管没有得到证实,为了获得具有稳定解的PML,可以使用等式中给出的 不等式。(8)或式中(9)应该满足,这意味着PML应用在任意各向异性介质 中的波模拟中受到限制。此外,值得注意的是,在违反第(7)种情况的情况下, 计算域中的各向异性波动方程将支持具有慢度矢量和群速度矢量相反方向的平 面波解,因此PML可以提供快速的指数增长解并变为在实践中不稳定。如 Bécache等人所示。甚至可以在入射波到达PML域的外部Dirichlet边界之前观 察到这种不稳定性。
接下来,本发明调查方程式的满意度。Meza-Fajardo&Papageorgiou(2008) 中总结的各向异性介质的(7),(8)或(9),为方便起见,在表1中也再次 列出。容易验证方程式。式(7)不适用于中等III,也不适用于磷灰石和锌,后 者的慢度和群速度图在图2中给出,表明慢度曲线不是凸的。都没有(8)或等 式(9)保持IV和V类型。仅对类型I,II和各向同性类型执行方程式。(7) 和等式(8)按住。因此,在将PML应用于材料III,磷灰石或锌时,可以观察 到快速的不稳定性,并且IV和V型的PML也不稳定。
表1
Figure BDA0002990702700000121
2、用于二阶线性弹性各向异性波系统的M-UPML或M-UFSPML
由于将PML解释为实际波动方程的解析延续复杂的副本的空间域,复杂的 坐标拉伸技术已经成为PML推导的主要工具。虽然PML是经典地通过以强形 式如式(1),也可以通过波浪的复数坐标拉伸得出用弱形式写的方程式。以经 典方式,在PML中波场方程的推导及其边界和/或界面条件是独立的,并且两者 可能不正确地匹配,这可能导致数值不稳定和数值精度的降低,如图3所示。 但是,采用弱形式的方法自然可以避免这种不匹配,因为自由和界面条件同时 以一致的方式拉伸。所获得的PML也是弱形式,因此可以通过FE/SE方法直 接离散化。在以下各节中,本发明将使用弱形式方法首先得出各向异性情况下 的二阶M-UPML或M-UFSPML。然后,本发明提供了一些公式与其他公式之间 的比较元素。
2.1弱形式M-UPML或M-UFSPML的推导
使用弱形式方法推导M-UPML或M-UFSPML的方法与强形式相同,只是 复数坐标拉伸函数的形式不同。对于M-UPML,使用经典的复数坐标拉伸函数, 而对于M-UFSPML,使用频移的函数。本发明从二阶形式开始形成弱形式的时 谐弹性动力学方程:
Figure BDA0002990702700000131
其中Ω表示PML域,Γ表示Ω在a上的边界符号表示其傅立叶变换。映 射方程(15)进入具有以下定义的复杂坐标的复杂空间:
Figure BDA0002990702700000132
Figure BDA0002990702700000133
本发明得到:
Figure BDA0002990702700000134
条件:
Figure BDA0002990702700000135
为了获得适合于实际空间中数值离散的PML公式,本发明需要对方程进行 变换,成实坐标(18)。根据衍生工具的链式规则和公式(16),(17),本 发明有:
Figure BDA0002990702700000136
Figure BDA0002990702700000137
代入式将(20)和(21)变成(18),本发明得到:
Figure BDA0002990702700000141
Figure BDA0002990702700000142
Figure BDA0002990702700000143
Figure BDA0002990702700000144
所得方程通过逆傅立叶变换回时域转换,可以写成:
Figure BDA0002990702700000145
可以在以下位置找到术语
Figure BDA0002990702700000146
的详细表达。很容易验证公式。给定的弱形式二阶PML由式(26)可以通过基于位移的数值方法直接离 散化。通过为
Figure BDA0002990702700000147
可以通过数值技术迭代地计算卷积总 结。等式(26)对于在只需将sx(x)或sy(y)设置为零,即可(仅在x或y方向 上)(即,在PML层的拐角之外)。
遵循Meza-Fajardo&Papageorgiou的想法,通过考虑类型转换,将提议的二 阶PML扩展为二阶M-PML:
Figure BDA0002990702700000148
Figure BDA0002990702700000149
以严格的方式,基于链式规则和等式(27),(28),本发明得到:
Figure BDA00029907027000001410
Figure BDA00029907027000001411
与:
Figure BDA0002990702700000151
(29)和等式(30)将导致M-PML的卷积积分数量是Meza-Fajardo& Papageorgiou(2008)中采用的近似链法则的两倍。与方程式相同(20)和(21), 只是将sx(x),sy(x)替换为sx(x,y),sy(x,y)。此外,在M-PML推导中采用严格 的链规则无法补偿完美匹配的特性的损失,已经通过在与该界面相切的方向上 引入额外的非零阻尼而恶化了。因此,为简化起见,本发明在推导中也使用了 近似链规则。通过弱形式波动方程的复拉伸,推导出2D和3D情况下的弱形式 二阶M-UFSPML。其派生的更多详细信息在附录B中给出。
2.2获得的二阶M-PML或M-UFSPML与其他配方的比较
新获得的M-UFSPML与得出的二阶UFSPML具有相同的表达。基于强形 式二阶波动方程得出的M-UFSPML或M-PML并非如此。在后一种情况下,通 常获得的二维情况下的频域二阶M-PML或M-UFSPML由下式给出:
Figure BDA0002990702700000152
Figure BDA0002990702700000153
Figure BDA0002990702700000154
Figure BDA0002990702700000155
sx(x,y)和sy(x,y)取决于x和y变量。因此,与二阶不同在各向同性的情况 下UFSPML推导,sx(x,y)×sy(x,y)的乘法不能同化
Figure BDA0002990702700000156
在 M-UFSPML中为
Figure BDA0002990702700000157
或M-UPML派生,使时域对应于M-UPML或M-UFSPML如果不引入复杂的基函数,就不能通过FE/SE直接离 散化。为了避免使用复杂的基函数,L提取所有的卷积项作为外力。例如,频移 坐标拉伸功能:
Figure BDA0002990702700000161
Figure BDA0002990702700000162
重写了(32)为:
Figure BDA0002990702700000163
Figure BDA0002990702700000164
然而,沿着不同媒介之间的界面,同化的分歧应力作为外力体项不能确保 正确计算
Figure BDA0002990702700000165
除了非常密集的网格外,这使其方法不适用于异质无界域 的截断。给出的一阶M-UFSPML存在相同的问题。
受到Komatitsch&Tromp想法的启发,通过变量拆分(例如拆分等式)获得 了准备用于FE/SE离散化的M-PML。公式(32):
Figure BDA0002990702700000166
Figure BDA0002990702700000167
Figure BDA0002990702700000168
=M-PML公式继承了Komatitsch&Tromp=的PML公式的一个缺点,由于位 移的分裂导致存在混合二阶和三阶微分方程,因此需要在内部区域修改数值算 法。时域对应物。此外,分裂位移变量的数量是原始变量的五倍=。
本发明的2D M-UFSPML涉及的最终方程式由14个方程式组成:2个方程 式用于位移向量,12个方程式用于存储变量。与经典的一阶PML或优化的一阶 拆分PML所需的20个方程相比,这在内存存储方面是一项重大改进。它与Martin,Komatitsch和Gedney给出的一阶M-UFSPML的公式相当,它由13个 方程组成:2个速度矢量方程,3个应力张量方程和8个记忆变量。对于二阶 M-PML,Ping等人的公式。由15个方程式组成:4个方程式,五个分裂位移矢量中的2个具有三阶时间导数,6个方程式,五个分裂位移矢量中的3个,具有 二阶时间导数。本发明的2D M-UPML由6个方程式组成:2个方程式用于位移 向量,4个方程式用于存储变量。
下面结合测试对本发明的技术效果作详细的描述。
1、现在以数字方式显示本发明的公式对于使用均质和分层非均质各向异性 介质截断无限域是稳定的。除非另有说明,否则在所有示例中,本发明对多项 式使用N=4的多项式每个光谱元素内的拉格朗日插值。根据Collino&Tsogka, 垂直于PML接口的阻尼曲线为:
Figure BDA0002990702700000171
其中L是PML的厚度,xi是从界面测量的沿法线方向的距离,而VP,max是 介质的最大P波速度。继Zhang&Shen之后,在M-UFSPML或UFSPML中, 本发明还设置了:
α(xi)=3πf0(1-|xi|/L), (44)
其中f0是源时间函数的主导频率:
Figure BDA0002990702700000172
其中t0是其开始时间。本发明研究整个域(包括PML域)中能量随时间的 衰减,以通过评估总能量(即动能和势能之和)来检查配方的长期稳定性:
Figure BDA0002990702700000173
其中
Figure BDA0002990702700000174
是速度向量的范数,而
Figure BDA0002990702700000175
是应力和应变组件。
1.1在各向同性介质的均匀半无限域中的波模拟
让本发明考虑在其顶部具有无牵引力边界的半无限均质各向同性域的截 断,而沿着其他三个侧面,该域被本发明的M-UFSPML截断。介质的弹性特性 和离散化参数为在表3的第二列中给出。对于M-UFSPML,本发明设置px|y和 py|xM-UFSPML的x等于0.1。为了比较吸收效率,四种无限域已经使用了截断 方法,即具有不同元素数量的M-UFSPML,M-UPML,Stacey吸收边界条件和 标准UFSPML吸收层,其优异的吸收效率已经得到了证明。除了α(xi)=0之外, M-PML的参数设置与M-UFSPML中的参数相同。在为了公平起见,在使用Stacey ABC的情况下,计算域将扩大以包括M-UFSPML占用的空间。在图3a, 3b,3c示出了由9层M-UFSPML,3层UFSPML和9层M-UPML在不同时间 获得的波场的快照。他们几乎一样。两者在数值上都是长期稳定的,如图2所 示。由图3e示出。3层,9层和15层M-UFSPML的吸收效率始终明显优于Stacey ABC。对于9层和15层M-PML,在吸收体波和瑞利波方面,吸收效率都接近3 层UFSPML。在使用相同数量的元素层的情况下,特别是对于近掠入射波, M-UFSPML的吸收效率明显好于M-UPML。
1.2均质各向异性无限域中的波模拟
现在,让本发明考虑使用本发明的M-UFSPML公式截断无限均质各向异性 域。各向异性介质的弹性特性和离散化参数在表2的第二和第三栏中给出。在 所有情况下,本发明设置px|y和py|xM-PML的x等于0.1,但III类各向异性介质 除外,对于本发明使用px|y=py|x=0.2。对于每种情况,本发明将9层M-UFSPML 的吸收效率与标准Stacey ABC的吸收效率进行比较。在图4至图8中,对于每 种情况,本发明首先显示由9层M-UFSPML获得的波场的快照,然后是总和动 能的时间演化曲线,最后归一化为水平和接收器上位移矢量的垂直分量。本发 明观察到相当不错9层M-UFSPML可以达到吸收效率,特别是对于材料V,用 于锌和磷灰石。但是,由于不完美,也可以观察到反射M-UFSPML的接口。
表2
Figure BDA0002990702700000191
1.3在各向同性和各向异性介质中分层均匀的波模拟。
在这一部分中,以展示M-UFSPML在异构介质中的性能。本发明考虑一个 二维的四层模型(图9),其介质属性在表3中给出。第一层和第四层的介质是 各向同性的。在第二层中,考虑VTI介质。为了考虑到地震学中使用的主要各 向异性介质,在第三层中将介质设置为各向同性介质中嵌入的一组垂直裂缝, 从而导致水平横向各向同性。的离散化参数由表4给出px|y和py|x如果 M-UFSPML的x等于0.1,则将M-UFSPML中包含的元素层数设置为9。
如图10所示,由于多次反射和透射,M-UFSPML在分层模型中吸收复杂波 的性能很好。可以实现长时间的数值稳定仿真,如图11所示,其物理域中的能 量衰减曲线说明了该问题。但是,从接收器100m靠近右侧区域的归一化位移的 两个分量的时间历史来看如图12所示,M-UFSPML可以观察到清晰的杂散反射 回到物理域。与使用Stacey ABC进行边界截断和扩展解所获得的结果相比, M-UFSPML的吸收效率比Stacey ABC更好。因此,本发明的M-UFSPML是在 实际应用中截断分层各向异性介质的非常有效的工具。
表3
Figure BDA0002990702700000201
表4
Figure BDA0002990702700000202
本发明的二阶未分裂M-UPML和M-UFSPML公式,用于笛卡尔域中的线性 弹性各向异性波的时域分析。由于本发明的公式是通过以弱形式编写的二阶各 向异性波系统的复杂拉伸得出的,因此它可以轻松地用于基于标准位移的有限 元或频谱元素波模拟。它易于实现,并允许自然应用无牵引力边界条件或界面 条件。此外,与其他现有公式相比,本发明的公式需要较少的变量数量和更少 的计算工作,这在解决大型三维问题时具有优势。通过使用具有二阶显式 Newmark时间积分器的高阶SE离散化,证明了该配方的使用,这表明了该方法 的准确性以及其解决正交各向异性对称各向异性介质中波传播的能力。尽管 M-PML的界面不是完美匹配,但与标准Stacey ABC等低阶吸收边界条件相比, 它在吸收体波和界面波方面仍然更加有效。未来研究的一个有趣的话题是将这 项工作扩展到各向异性粘弹性介质中的前向波和伴随波模拟。
以上所述仅为本发明的较佳实施例而已,并不用以限制本发明,凡在本发 明的精神和原则之内所作的任何修改、等同替换和改进等,均应包含在本发明 的保护范围之内。

Claims (10)

1.一种基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法包括以下步骤:
步骤一,对时域二阶线弹性各向异性波动方程的弱形式进行傅里叶变换,得到对应的频域方程;
步骤二,在多轴完全匹配层的概念基础上,采用复坐标延拓技术,对弱形式频域波动方程中的空间变量进行延拓,即得新构建的人工边界。
2.如权利要求1所述的基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述步骤一中多项式相乘的傅里叶反变换的方法首先通过微调分母参数以规避极点,并将相乘化为相加的形式,直接通过傅里叶反变换得到其时域表达式。
3.如权利要求1所述的基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述步骤二中复坐标延拓函数选用复频移形式。
4.如权利要求1所述的基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述步骤二中新构建的人工边界被称为非分裂卷积多轴完美匹配层。
5.如权利要求1所述的基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法还包括以下步骤:
第一步,通过对弱形式的基于位移的各向异性波方程进行复杂的拉伸,得到新型的时域二阶多轴未分裂PML或多轴未分裂频移PML使用经典或频移复数坐标拉伸函数得出;
第二步,获得的时域二阶多轴未分裂PML或多轴未分裂频移PML以位移为基础的弱形式,使用有限元或频谱元素方法进行少量修改实现。
6.如权利要求5所述的基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法的异性波动方程为:遵循使用的矩阵表示法,各向异性介质由动量方程控制:
Figure FDA0002990702690000021
和应力-应变关系:
σ=C·e (2)
其中u=[ux uy uz]T=[u1 u2 u3]T和f=[f1 f2 f3]T是位移,身体力量向量:
σ=[σxx σyy σzz σyz σxz σxy]T=[σ1 σ2 σ3 σ4 σ5 σ6]T
e=[exx eyy ezz exx exz exy]T=[e1 e2 e3 e4 e5 e6]T
Figure FDA0002990702690000022
在最一般的各向异性情况下,具有三斜对称性的弹性模型,弹性矩阵为:
Figure FDA0002990702690000023
cmw(m,n=1,2,…,6)是由四阶张量Ciijk(i,j,k,l=1,2,3)形成的按照以下方式:m对应于第一对索引i,j,n对应于第二对索引k,1,对应m→i,j和n→k,l如下:1→1,1;2→2,2;3→3,3;4→2,3;5→1,3;6→1,2;垂直横向各向同性VTI对称性,水平横向各向同性HTI对称和六角对称各向异性模型的弹性矩阵正交各向异性对称模型读取:
Figure FDA0002990702690000031
在二维x-y笛卡尔坐标情况下,通过以下方法获得面内波动方程,有u=[u1 u2]T,σ=[σ1σ2 σ6]T,e=[e1 e2 e6]T
Figure FDA0002990702690000032
7.如权利要求5所述的基于地震各向异性二阶波系统的M-UFSPML模型构建方法,其特征在于,所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法的通过对PML接受的平面波解的扰动分析,分别获得了必要的和充分的稳定性条件,为了与(6)中给出的弹性矩阵表示法一致,在x方向上吸收的PML的必要稳定性条件为:
[(c12+c66)2-c11(c22-c66)]×[(c12+c66)2+c66(c22-c66)]≤0 (7)
足够的稳定性条件表示,如果实现以下两个条件之一,则足以满足PML系统在x方向上的吸收的稳定性:
(c12+c66)2≤(c11-c66)(c22-c66) (8)
Figure FDA0002990702690000033
8.一种由权利要求1~7任意一项所述基于地震各向异性二阶波系统的M-UFSPML模型构建方法构建的M-UFSPML模型。
9.一种搭载权利要求8所述M-UFSPML模型的智能终端。
10.一种搭载权利要求8所述M-UFSPML模型的服务器。
CN202110314751.0A 2021-03-24 2021-03-24 M-ufspml模型、构建方法、智能终端、服务器 Active CN113435074B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110314751.0A CN113435074B (zh) 2021-03-24 2021-03-24 M-ufspml模型、构建方法、智能终端、服务器

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110314751.0A CN113435074B (zh) 2021-03-24 2021-03-24 M-ufspml模型、构建方法、智能终端、服务器

Publications (2)

Publication Number Publication Date
CN113435074A true CN113435074A (zh) 2021-09-24
CN113435074B CN113435074B (zh) 2024-02-09

Family

ID=77752902

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110314751.0A Active CN113435074B (zh) 2021-03-24 2021-03-24 M-ufspml模型、构建方法、智能终端、服务器

Country Status (1)

Country Link
CN (1) CN113435074B (zh)

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102880590A (zh) * 2012-09-25 2013-01-16 复旦大学 二阶波动方程的非分裂完全匹配层的构造方法
US20150127311A1 (en) * 2013-11-06 2015-05-07 Weidlinger Associates, Inc. Computer Implemented Apparatus and Method for Finite Element Modeling Using Hybrid Absorbing Element
CN106547023A (zh) * 2017-01-16 2017-03-29 青岛海洋地质研究所 一种精度高、计算稳定的复杂介质地震波场延拓方法
CN109188512A (zh) * 2018-09-17 2019-01-11 中国石油大学(华东) 基于非规则扇形网格剖分的起伏隧道空间正演模拟系统及方法
CN110058303A (zh) * 2019-05-06 2019-07-26 吉林大学 声波各向异性逆时偏移混合方法
WO2020063131A1 (zh) * 2018-09-28 2020-04-02 中国科学院地质与地球物理研究所 拟微分算子储存方法及装置

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102880590A (zh) * 2012-09-25 2013-01-16 复旦大学 二阶波动方程的非分裂完全匹配层的构造方法
US20150127311A1 (en) * 2013-11-06 2015-05-07 Weidlinger Associates, Inc. Computer Implemented Apparatus and Method for Finite Element Modeling Using Hybrid Absorbing Element
CN106547023A (zh) * 2017-01-16 2017-03-29 青岛海洋地质研究所 一种精度高、计算稳定的复杂介质地震波场延拓方法
CN109188512A (zh) * 2018-09-17 2019-01-11 中国石油大学(华东) 基于非规则扇形网格剖分的起伏隧道空间正演模拟系统及方法
WO2020063131A1 (zh) * 2018-09-28 2020-04-02 中国科学院地质与地球物理研究所 拟微分算子储存方法及装置
CN110058303A (zh) * 2019-05-06 2019-07-26 吉林大学 声波各向异性逆时偏移混合方法

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
ZHINAN XIE, DIMITRI KOMATITSCH, ROLAND MARTIN AND REN´E MATZEN: "Improved forward wave propagation and adjoint-based sensitivity kernel calculations using a numerically stable finite-element PML", GEOPHYSICAL JOURNAL INTERNATIONAL, vol. 198, no. 3 *
ZHINAN XIE; RENÉ MATZEN; PAUL CRISTINI; DIMITRI KOMATITSCH; ROLAND MARTIN: "A perfectly matched layer for fluid-solid problems: Application to ocean-acoustics simulations with solid ocean bottoms", J. ACOUST. SOC. AM., vol. 140, no. 1, XP012209139, DOI: 10.1121/1.4954736 *
章旭斌;谢志南;廖振鹏: "SH波动模拟中透射边界反射放大失稳研究", 哈尔滨工程大学学报, vol. 40, no. 6 *
谢志南;章旭斌;: "弱形式时域完美匹配层", 地球物理学报, no. 10 *
谢志南;郑永路;章旭斌;唐丽华: "弱形式时域完美匹配层——滞弹性近场波动数值模拟", 地球物理学报, vol. 62, no. 8 *

Also Published As

Publication number Publication date
CN113435074B (zh) 2024-02-09

Similar Documents

Publication Publication Date Title
Carcione Wave fields in real media: Wave propagation in anisotropic, anelastic, porous and electromagnetic media
Zhang et al. Unsplit complex frequency-shifted PML implementation using auxiliary differential equations for seismic wave modeling
Gosselin-Cliche et al. 3D frequency-domain finite-difference viscoelastic-wave modeling using weighted average 27-point operators with optimal coefficients
Askan et al. Full waveform inversion for seismic velocity and anelastic losses in heterogeneous structures
Alkhalifah Scanning anisotropy parameters in complex media
Liu et al. A mixed-grid finite element method with PML absorbing boundary conditions for seismic wave modelling
Marburg Boundary element method for time-harmonic acoustic problems
Guzina et al. On the analysis of wave motions in a multi-layered solid
Liu et al. A new kind of optimal second-order symplectic scheme for seismic wave simulations
Hao et al. Viscoacoustic anisotropic wave equations
CN109946742A (zh) 一种TTI介质中纯qP波地震数据模拟方法
Zhuang et al. A simple implementation of PML for second-order elastic wave equations
Yang et al. An optimal nearly analytic discrete-weighted Runge-Kutta discontinuous Galerkin hybrid method for acoustic wavefield modeling
Day et al. Seismic response of hemispherical foundation
Xu et al. Time-space-domain temporal high-order staggered-grid finite-difference schemes by combining orthogonality and pyramid stencils for 3D elastic-wave propagation
CN109188517B (zh) 基于Higdon余弦型加权的混合吸收边界条件方法
CN109725346B (zh) 一种时间-空间高阶vti矩形网格有限差分方法及装置
CN113435074A (zh) M-ufspml模型、构建方法、智能终端、服务器
Brick et al. Fast direct solution of 3-D scattering problems via nonuniform grid-based matrix compression
He et al. Modeling 3-D elastic wave propagation in TI media using discontinuous Galerkin method on tetrahedral meshes
Zhang et al. Viscoelastic wave propagation in transversely isotropic media based on constant-order fractional polynomial approximations
AlSalem et al. Embedded boundary methods for modeling 3D finite-difference Laplace-Fourier domain acoustic-wave equation with free-surface topography
Wang et al. Improved hybrid iterative optimization method for seismic full waveform inversion
Singla et al. On planar seismic wavefront modeling for estimating rotational ground motions: Case of 2-D SH line-source
Fu et al. An efficient high-order multiscale finite element method for frequency-domain elastic wave modeling

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