CN112163332A - 基于自然单元无限元耦合的2.5维自然电场数值模拟方法 - Google Patents

基于自然单元无限元耦合的2.5维自然电场数值模拟方法 Download PDF

Info

Publication number
CN112163332A
CN112163332A CN202011017516.9A CN202011017516A CN112163332A CN 112163332 A CN112163332 A CN 112163332A CN 202011017516 A CN202011017516 A CN 202011017516A CN 112163332 A CN112163332 A CN 112163332A
Authority
CN
China
Prior art keywords
natural
unit
wave number
infinite
dimensional
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.)
Pending
Application number
CN202011017516.9A
Other languages
English (en)
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.)
Central South University
Original Assignee
Central South University
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 Central South University filed Critical Central South University
Priority to CN202011017516.9A priority Critical patent/CN112163332A/zh
Publication of CN112163332A publication Critical patent/CN112163332A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/15—Correlation function computation including computation of convolution operations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06Q—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES; SYSTEMS OR METHODS SPECIALLY ADAPTED FOR ADMINISTRATIVE, COMMERCIAL, FINANCIAL, MANAGERIAL OR SUPERVISORY PURPOSES, NOT OTHERWISE PROVIDED FOR
    • G06Q50/00—Information and communication technology [ICT] specially adapted for implementation of business processes of specific business sectors, e.g. utilities or tourism
    • G06Q50/06—Energy or water supply
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/04—Constraint-based CAD
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Business, Economics & Management (AREA)
  • Health & Medical Sciences (AREA)
  • Economics (AREA)
  • Evolutionary Computation (AREA)
  • Operations Research (AREA)
  • Geometry (AREA)
  • Computer Hardware Design (AREA)
  • Computing Systems (AREA)
  • Public Health (AREA)
  • Water Supply & Treatment (AREA)
  • General Health & Medical Sciences (AREA)
  • Human Resources & Organizations (AREA)
  • Marketing (AREA)
  • Primary Health Care (AREA)
  • Strategic Management (AREA)
  • Tourism & Hospitality (AREA)
  • General Business, Economics & Management (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本发明提供了一种基于自然单元无限元耦合的2.5维自然电场数值模拟方法，具体步骤如下：S1、建立二维自然电场地电模型，依据模型参数在自然单元区域设置离散点集，在区域边界扩展无限单元充当边界单元；S2、构建2.5维自然电场地电模型满足的微分方程，推导其基于自然单元无限元耦合的变分形式及其积分形式；S3、构建基于Laplace插值函数的自然单元法，推导2.5维自然单元形函数及其导数的基本方程；S4、构建2.5维无限元形函数及其导数的基本方程；S5、计算各波数条件下的总刚度矩阵；S6、处理电源信息；S7、求解大型稀疏方程组得到各波数条件下的波数域自然电位剖面，最后将各波数域自然电位剖面通过傅里叶反变换得到空间域下的自然电位剖面。

Description

基于自然单元无限元耦合的2.5维自然电场数值模拟方法
技术领域
本发明涉及一种基于自然单元无限元耦合的2.5维自然电场数值模拟算法，特别是涉及一种适用于工程与环境相关的复杂、多源、动态源模型的自然电场高精度高效率正演模拟的方法。
背景技术
自然电场法作为一种被动源法，对多种与污染扩散伴生的渗流运动、氧化还原以及微生物活动信号较为敏感，适用于工程与环境物探领域，尤其是地下有机污染物检测及监测。
当前应用于稳定电流场数值模拟的方法各有优势，但也存在一定的缺陷。有网格的数值方法如有限元、有限差分等，需要构建网格单元，该过程不仅增加了前期处理工作量，也使得复杂模型数值模拟变得困难；而无单元法虽不受网格约束，但计算量巨大，且其形函数不满足插值性质，造成本质边界条件施加困难。
自然单元法依托于求解域内离散节点的Voronoi结构，采用无网格方式构造插值函数，免去了有限元法需手动或半自动剖分网格的麻烦，且模型不受网格约束，适合模拟复杂地电模型。同时针对在处理地电模型中无限域或半无限域问题时所引起的截断边界误差，引入无限单元，实现其在区域边界与自然单元节点的有效耦合，形成一种适用于稳定电流场2.5维高精度高效率数值模拟的新方法。
发明内容
本发明目的在于提供一种基于Laplace插值算法，构造自然单元形函数，同时引入无限元作为边界单元，解决无限元或半无限域所引起的截断边界误差，实现基于自然单元无限元耦合的2.5维自然电场数值模拟的方法。
为实现上述目的，本发明提供了一种基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其具体步骤如下：
S1：建立二维自然电场地电模型，依据模型参数在自然单元区域设置离散点集，在区域边界扩展无限单元充当边界单元；
S2：构建2.5维自然电场地电模型满足的微分方程，推导其基于自然单元无限元耦合的变分形式以及积分形式；
S3：构建基于Laplace插值函数的自然单元法，推导二维自然单元形函数及其导数的基本方程；
S4：构建二维无限元形函数及其导数的基本方程；
S5：基于S2、S3、S4推导的公式，计算各波数条件下的总刚度矩阵；
S6：处理电源信息；
S7：求解大型稀疏方程组得到各波数条件下的波数域自然电位剖面，然后将各波数域自然电位剖面通过傅里叶反变换得到空间域下的自然电位剖面。
作为本发明的进一步方案：所述步骤S1中所建立的二维自然电场实施例地电模型，其模型中自然单元区域尺度设置为20m×30m(深度×横向)，布设651个节点。
所述二维自然电场地电模型除地表外的三个边界由无限单元作为边界单元充当边界条件，其总单元数为70；设定地下水位位于2m深度处，上层电阻率为20欧姆米，下层电阻率为100欧姆米；通过离散正负点电源近似正负电荷的区域性分布，设定各点电源幅值为2mA。
所述自然单元区域由Laplace自然邻点插值函数求取，边界区域由无限单元形函数求取。
作为本发明的进一步方案：所述步骤S3中的自然单元法采用的是一种基于Voronoi图和Delaunay三角化几何结构，以自然邻点插值为试函数的一种无网格数值方法；
作为本发明的进一步方案：所述步骤S5中计算各波数条件下的总刚度矩阵的具体过程如下：
(1)对自然单元区域的各背景积分网格进行循环；
(2)对各自然单元区域的背景积分网格中的高斯积分点进行循环；
(3)搜索自然单元区域的背景积分网格中的各高斯积分点的自然邻点；
(4)计算自然单元区域的背景积分网格中的各自然邻点相关于高斯积分点的形函数及导数；
(5)计算自然单元区域的背景积分网格中的各高斯积分点的子刚度矩阵；
(6)将自然单元区域的背景积分网格中的各子刚度矩阵加载到总刚度矩阵中；
(7)同时在各波数循环条件下，对边界无限单元进行循环；
(8)计算各无限单元子刚度矩阵；
(9)将各无限单元子刚度矩阵加载到总刚度矩阵中。
所述波数及其权值选择如下：
k＝[0.004758,0.0407011,0.1408855,0.393225,1.088038]；
g＝[0.0099472,0.031619,0.0980327,0.2511531,0.7260814。.
作为本发明的进一步方案：所述步骤S7中的具体计算过程为：
(1)在各波数条件下，求解得到当前波数对应的总刚度矩阵；
(2)随后联合电源项，求解大型稀疏方程组，得到各波数条件下的波数域总自然电位剖面；
(3)最后通过傅里叶反变换计算得到空间域中的总自然电位剖面。
应用本发明的技术方案，具有以下有益效果：
(1)本发明提供的方法，实现基于Laplace自然邻点插值的自然单元与无限单元在空间及物性上的有效耦合，并通过傅里叶正反变换，实现波数域到空间域的2.5维自然电场数值模拟。
(2)本发明提供的方法，不受网格单元的约束，节点布置灵活，适用于物性复杂的模型。
(3)本发明提供的方法，在模型边界处耦合无限单元，能有效解决截断边界问题，模拟精度高，计算效率高。
(4)本发明提供的方法，自然单元区域的计算最终都可转化为简单的代数运算，计算速度快，编程实现简便。
(5)本发明提供的方法，能推动自然电场等多源、动态源复杂地电模型数值模拟工作的开展，为相应地球物理方法应用到污染监测、检测以及其他工程与环境地球物理问题中提供数值计算基础。
(6)本发明中通过采用基于Voronoi图和Delaunay三角化几何结构、以自然邻点插值为试函数的一种无网格数值方法推导2.5维自然单元形函数及其导数的基本方程，使其具有无网格法和有限元法的优点，又克服了两者的缺陷。
(7)本发明中对电源信息的处理，采用离散点源近似模拟呈区域分布的自然电场源，需在相应节点位置赋点源幅值；同时因无限元的引入，故不需考虑源分布对边界条件的影响，适用于多源、动态源模型。
(8)本发明通过计算方式将各波数域自然电位剖面通过傅里叶反变换得到空间域下的自然电位剖面，其数值模拟结果表明，地质微生物降解地下有机污染物时所发生的氧化还原反应能在地表产生负的自然电位异常；同时表明本发明所提的数值模拟方法是有效的可行性。
除了上面所描述的目的、特征和优点之外，本发明还有其它的目的、特征和优点。下面将参照图，对本发明作进一步详细的说明。
附图说明
构成本申请的一部分的附图用来提供对本发明的进一步理解，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：
图1是本发明基于自然单元无限元耦合的2.5维自然电场地电模型示意图；
图2是本发明中节点A的Voronoi单胞示意图；
图3是本发明中7节点Voronoi图；
图4是本发明中Delaunay外接圆的示意图；
图5是本发明中Delaunay三角化的示意图；
图6是本发明中7节点Voronoi图中的x计算点的示意图；
图7是本发明中x计算点的自然邻点及其Laplace插值元素的示意图；
图8是本发明中二维四节点无限元映射过程示意图；
图9是本发明中实施例地表自然电位异常示意图。
具体实施方式
以下结合附图对本发明的实施例进行详细说明，但是本发明可以根据权利要求限定和覆盖的多种不同方式实施。
参见图1至图9所示，本发明提供的一种基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其具体步骤如下：
S1：建立二维自然电场地电模型，依据模型参数在自然单元区域设置离散点集，在区域边界扩展无限单元充当边界单元；
S2：构建2.5维自然电场地电模型满足的微分方程，推导其基于自然单元无限元耦合的变分形式以及积分形式；
S3：构建基于Laplace插值函数的自然单元法，推导二维自然单元形函数及其导数的基本方程；
S4：构建二维无限元形函数及其导数的基本方程；
S5：基于S2推导的变分问题，通过S3所推导的公式可计算得到K1，通过S4所推导的公式可计算得到K2，以计算各波数条件下的总刚度矩阵；
S6：处理电源信息；
S7：求解大型稀疏方程组得到各波数条件下的波数域自然电位剖面，然后将各波数域自然电位剖面通过傅里叶反变换得到空间域下的自然电位剖面。
作为本发明的进一步方案：所述步骤S1中所建立的二维自然电场地电模型，其模型中自然单元区域尺度设置为20m×30m(深度×横向)，布设651个节点。
所述二维自然电场地电模型除地表外的三个边界由无限单元作为边界单元充当边界条件，其总单元数为70。
设定地下水位位于2m深度处，上层电阻率为20欧姆米，下层电阻率为100欧姆米。
通过离散正负点电源近似正负电荷的区域性分布，设定各点电源幅值为2mA。
作为本发明的进一步方案：所述步骤S2中2.5维自然电场波数域总电位U的边值问题的表达示为：
Figure BDA0002699564650000051
Figure BDA0002699564650000052
其中：σ为介质电导率，I0为点源电流强度，δ(A)为与点源位置A有关的δ函数，n为边界外法向单位矢量，
Figure BDA0002699564650000053
为2维哈密顿算子、且
Figure BDA0002699564650000054
k为离散波数，U为三维空间中的电位u沿y方向进行傅里叶变换的波数域总电位，Ω、Γs分别为经傅里叶变换后的研究区域、地表边界。
则，所述2.5维自然电场波数域总电位U的边值问题所对应的变分问题的表达示为：
Figure BDA0002699564650000055
δF(U)＝0 (4)
将式(3)右端进行变分，则为：
Figure BDA0002699564650000056
假设自然单元区域的某一计算点包含n个自然邻点，利用该n个节点构造自然单元形函数，采用该组形函数φi作为试函数近似计算点x＝(x,z)电位U：
Figure BDA0002699564650000057
将(6)式代入(5)式有：
Figure BDA0002699564650000061
对(7)式中的面积分项做变换有：
Figure BDA0002699564650000062
所述无限元区域所满足的偏微分方程及其变分形式，同样可分别由式(1)、(2)和(8)表示。故对于无限元区域的面积分的表达示为：
Figure BDA0002699564650000063
其中：式(8)和式(9)中的Φ分别由自然单元及无限元求取，合并式(8)和式(9)，有：
δUT[(K1+K2)U-F]＝0 (10)
其中有：
Figure BDA0002699564650000064
Figure BDA0002699564650000065
F＝∫ΩI0δ(A)ΦdΩ (13)
由于δU是任意的，故(10)式成立的条件为：
(K1+K2)U＝F (14)
即：
KU＝F (15)
其中：K＝K1+K2为总系数矩阵，而自然单元和无限元耦合的具体过程可表示为：
Figure BDA0002699564650000066
其中：Kij为耦合方法的总刚度矩阵元素，
Figure BDA0002699564650000071
为自然单元刚度元素，
Figure BDA0002699564650000072
为无限元刚度元素，i、j为节点编号。
无限元区域以无限单元作为积分单元计算K2；而自然单元区域可用一组背景单元进行剖分(如矩形网格、Delaunay三角网)，在背景单元内采用高斯积分计算K1。对于单个积分单元Ωe可写为：
Figure BDA0002699564650000073
Figure BDA0002699564650000074
Fe分别为自然单元积分单元Ωe的面积分子系数矩阵、无限元积分单元Ωe的面积分子系数矩阵、子右端项，可表示为：
Figure BDA0002699564650000075
Figure BDA0002699564650000076
Figure BDA0002699564650000077
假设背景单元Ωe中采用Ng个高斯积分点xg＝(xg,zg)，g＝1,2,…,Ng进行积分计算，它们对应的权重wg，g＝1,2,…,Ng和雅克比值Jg，g＝1,2,…,Ng，高斯点支持域内包含n个自然邻点，则有：
Figure BDA0002699564650000078
其中：
Figure BDA0002699564650000079
B＝σk2Φ(xg)ΦT(xg) (23)
展开(22)、(23)有：
Figure BDA00026995646500000710
Figure BDA0002699564650000081
式(24)和式(25)中的形函数及其导数，在不同区域用不同算法求取，即自然单元区域由Laplace自然邻点插值函数求取，边界区域由无限单元形函数求取。
作为本发明的进一步方案：所述步骤S3中的自然单元法采用的是一种基于Voronoi图和Delaunay三角化几何结构，以自然邻点插值为试函数的一种无网格数值方法；
对于二维空间，Laplace插值函数的形式为：
Figure BDA0002699564650000082
其中，
Figure BDA0002699564650000083
si(x,y)是与节点i关联的Voronoi边的长度，hi(x,y)是插值点到节点i的Voronoi边的垂直距离，si(x,y)、hi(x,y)均为节点i的坐标(x,y)的函数。
由Voronoi结构的定义可知，待求点的Voronoi结构顶点是由待求点x与相邻两自然邻节点组成的三角形外接圆圆心，si(x,y)实际是圆心间的距离；对于任意三角形ΔABC，设其顶点坐标为A(x1,y1)、B(x2,y2)、C(x3,y3)，则该三角形外接圆圆心V＝(vx,vy)可表示为(设C为计算点，即高斯积分点)：
Figure BDA0002699564650000084
Figure BDA0002699564650000085
Figure BDA0002699564650000086
Figure BDA0002699564650000087
Figure BDA0002699564650000088
Figure BDA0002699564650000091
其中
D＝2[(x1-x3)(y2-y3)-(x2-x3)(y1-y3)] (33)
α＝(x2+x1)(x2-x1)+(y2+y1)(y2-y1) (34)
D,x＝2(y1-y2) (35)
D,y＝2(x2-x1) (36)设任意一条Voronoi边的两个端点为V1＝(vx1,vy1)，V2＝(vx2,vy2)，则有
Figure BDA0002699564650000092
由式(26)可得φi(x,y)关于点i坐标的一阶导数为：
Figure BDA0002699564650000093
Figure BDA0002699564650000094
其中αi关于坐标(x,y)的一阶导数为：
Figure BDA0002699564650000095
Figure BDA0002699564650000096
其中有：
Figure BDA0002699564650000097
Figure BDA0002699564650000098
Figure BDA0002699564650000101
Figure BDA0002699564650000102
Figure BDA0002699564650000103
作为本发明的进一步实施例：在图1地电模型中，已在模型边界处注明无限单元的分布特征。本发明采用的是Astley型映射无限单元，其通过坐标映射使单元沿某一方向无限延伸，具体映射过程如图8所示。在二维四节点无限单元中，任意位置的空间坐标可表示为：
Figure BDA0002699564650000104
Figure BDA0002699564650000105
其中Mi为坐标映射函数，其表达式为：
Figure BDA0002699564650000106
Figure BDA0002699564650000107
Figure BDA0002699564650000108
Figure BDA0002699564650000109
式中ξi，ηi为局部坐标。经过上述坐标映射之后，电位借助无限单元延伸至无限远处，并且在无限远处衰减为零。
无限单元中任意位置的波数域总电位可表示为：
Figure BDA00026995646500001010
其中Ui为波数域中各节点电位；Ni为无限元插值形函数，其表达式为：
Figure BDA0002699564650000111
Figure BDA0002699564650000112
Figure BDA0002699564650000113
Figure BDA0002699564650000114
在无限单元中，雅克比矩阵[J]由映射函数求得，具体表达式为：
Figure BDA0002699564650000115
进一步地，无限元形函数对坐标的导数可表示为：
Figure BDA0002699564650000116
作为本发明的进一步方案：所述步骤S5中计算各波数条件下的总刚度矩阵的具体过程如下：
(1)对自然单元区域的各背景积分网格进行循环；
(2)对各自然单元区域的背景积分网格中的高斯积分点进行循环；
(3)搜索自然单元区域的背景积分网格中的各高斯积分点的自然邻点；
(4)计算自然单元区域的背景积分网格中的各自然邻点相关于高斯积分点的形函数及导数；
(5)计算自然单元区域的背景积分网格中的各高斯积分点的子刚度矩阵；
(6)将自然单元区域的背景积分网格中的各子刚度矩阵加载到总刚度矩阵中；
(7)同时在各波数循环条件下，对边界无限单元进行循环；
(8)计算各无限单元子刚度矩阵；
(9)将各无限单元子刚度矩阵加载到总刚度矩阵中。
优选的，所述波数及其权值选择如下：
波数k＝[0.004758,0.0407011,0.1408855,0.393225,1.088038]；
波数权值g＝[0.0099472,0.031619,0.0980327,0.2511531,0.7260814。.
作为本发明的进一步方案：所述步骤S6中的电源信息的电源项如式(20)所示，且其当形函数具有Kronecker delta函数性质时，源项积分为：
Figure BDA0002699564650000121
作为本发明的进一步方案：所述步骤S7中的具体计算过程为：
(1)在各波数条件下，求解得到当前波数对应的总刚度矩阵；
(2)随后联合电源项，求解大型稀疏方程组，得到各波数条件下的波数域总自然电位剖面；
(3)最后通过傅里叶反变换计算得到空间域中的总自然电位剖面。
以上所述仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (8)

1.一种基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于，具体步骤如下：
S1：建立二维自然电场地电模型，依据模型参数在自然单元区域设置离散点集，在区域边界扩展无限单元充当边界单元；
S2：构建2.5维自然电场地电模型满足的微分方程，推导其基于自然单元无限元耦合的变分形式以及积分形式；
S3：构建基于Laplace插值函数的自然单元法，推导二维自然单元形函数及其导数的基本方程；
S4：构建二维无限元形函数及其导数的基本方程；
S5：基于S2、S3、S4推导的公式，计算各波数条件下的总刚度矩阵；
S6：处理电源信息；
S7：求解大型稀疏方程组得到各波数条件下的波数域自然电位剖面，然后将各波数域自然电位剖面通过傅里叶反变换得到空间域下的自然电位剖面。
2.根据权利要求1所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：所述步骤S1中所建立的二维自然电场地电模型，其模型中自然单元区域的深度尺度设置为20m、横向尺度设置为30m，布设651个节点。
3.根据权利要求2所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：通过离散正负点电源近似地下正负电荷的区域性分布，设定各点电源幅值为2mA。
4.根据权利要求3所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：所述自然单元区域由Laplace自然邻点插值函数求取，边界区域由无限单元形函数求取。
5.根据权利要求1所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：所述步骤S3中的自然单元法采用的是一种基于Voronoi图和Delaunay三角化几何结构，以自然邻点插值为试函数的一种无网格数值方法。
6.根据权利要求1所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：所述步骤S5中计算各波数条件下的总刚度矩阵的具体过程如下：
(1)对自然单元区域的各背景积分网格进行循环；
(2)对各自然单元区域的背景积分网格中的高斯积分点进行循环；
(3)搜索自然单元区域的背景积分网格中的各高斯积分点的自然邻点；
(4)计算自然单元区域的背景积分网格中的各自然邻点相关于高斯积分点的形函数及导数；
(5)计算自然单元区域的背景积分网格中的各高斯积分点的子刚度矩阵；
(6)将自然单元区域的背景积分网格中的各子刚度矩阵加载到总刚度矩阵中；
(7)同时在各波数循环条件下，对边界无限单元进行循环；
(8)计算各无限单元子刚度矩阵；
(9)将各无限单元子刚度矩阵加载到总刚度矩阵中。
7.根据权利要求6所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：所述波数及其权值选择如下：
k＝[0.004758,0.0407011,0.1408855,0.393225,1.088038]；
g＝[0.0099472,0.031619,0.0980327,0.2511531,0.7260814]；
其中：k为波数，g为波数权值。
8.根据权利要求1所述的基于自然单元无限元耦合的2.5维自然电场数值模拟方法，其特征在于：所述步骤S7中的具体计算过程为：
(1)在各波数条件下，求解得到当前波数对应的总刚度矩阵；
(2)随后联合电源项，求解大型稀疏方程组，得到各波数条件下的波数域总自然电位剖面；
(3)最后通过傅里叶反变换计算得到空间域中的总自然电位剖面。
CN202011017516.9A 2020-09-24 2020-09-24 基于自然单元无限元耦合的2.5维自然电场数值模拟方法 Pending CN112163332A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202011017516.9A CN112163332A (zh) 2020-09-24 2020-09-24 基于自然单元无限元耦合的2.5维自然电场数值模拟方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202011017516.9A CN112163332A (zh) 2020-09-24 2020-09-24 基于自然单元无限元耦合的2.5维自然电场数值模拟方法

Publications (1)

Publication Number Publication Date
CN112163332A true CN112163332A (zh) 2021-01-01

Family

ID=73863841

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202011017516.9A Pending CN112163332A (zh) 2020-09-24 2020-09-24 基于自然单元无限元耦合的2.5维自然电场数值模拟方法

Country Status (1)

Country Link
CN (1) CN112163332A (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112163354A (zh) * 2020-09-24 2021-01-01 中南大学 一种基于自然单元的2.5维自然电场数值模拟方法
CN113221409A (zh) * 2021-05-07 2021-08-06 桂林理工大学 一种有限元和边界元耦合的声波二维数值模拟方法和装置
CN116227308A (zh) * 2023-05-09 2023-06-06 广东石油化工学院 一种浅层测井自然电场的数值模拟方法及系统

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109740230A (zh) * 2018-12-26 2019-05-10 中南大学 一种自然电场三维多向映射耦合数值模拟方法

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109740230A (zh) * 2018-12-26 2019-05-10 中南大学 一种自然电场三维多向映射耦合数值模拟方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
JING XIE 等: "2.5D self-potential forward modeling by natural-infinite element coupling method", 《JOURNAL OF APPLIED GEOPHYSICS》 *
麻昌英 等: "基于全局弱式无单元法直流电阻率正演模拟", 《地球物理学报》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112163354A (zh) * 2020-09-24 2021-01-01 中南大学 一种基于自然单元的2.5维自然电场数值模拟方法
CN113221409A (zh) * 2021-05-07 2021-08-06 桂林理工大学 一种有限元和边界元耦合的声波二维数值模拟方法和装置
CN116227308A (zh) * 2023-05-09 2023-06-06 广东石油化工学院 一种浅层测井自然电场的数值模拟方法及系统

Similar Documents

Publication Publication Date Title
CN112163332A (zh) 基于自然单元无限元耦合的2.5维自然电场数值模拟方法
CN106199742B (zh) 一种频率域航空电磁法2.5维带地形反演方法
Chen et al. A method of DEM construction and related error analysis
Spitzer A 3-D finite-difference algorithm for DC resistivity modelling using conjugate gradient methods
CN105426339B (zh) 一种基于无网格法的线源时域电磁响应数值计算方法
Montero et al. A 3-D diagnostic model for wind field adjustment
CN106980736B (zh) 一种各向异性介质的海洋可控源电磁法有限元正演方法
CN103728667B (zh) 一种视三维高密度电法的浅表层地质结构建模方法
CN108388707B (zh) 一种三维非对称结构土壤模型下基于场路耦合的直流偏磁计算方法
CN108108579B (zh) 直流电阻率无单元法中耦合有限单元法的边界处理方法
CN113553748B (zh) 一种三维大地电磁正演数值模拟方法
CN111638556B (zh) 基于地空分解策略的大地电磁正演方法及装置、存储介质
CN112949134A (zh) 基于非结构有限元方法的地-井瞬变电磁反演方法
CN115238550A (zh) 自适应非结构网格的滑坡降雨的地电场数值模拟计算方法
CN112307640A (zh) 基于自然单元-无限元的三维多源自然电位数值模拟方法
Tossavainen et al. Estimating shapes and free surfaces with electrical impedance tomography
CN112163354A (zh) 一种基于自然单元的2.5维自然电场数值模拟方法
CN109101463A (zh) 一种广域层状大地格林函数的多精度求解方法
Akça et al. Extraction of structure-based geoelectric models by hybrid genetic algorithms
CN103698810A (zh) 混合网最小走时射线追踪层析成像方法
Zhang et al. Wind modelling for wind erosion research by open source computational fluid dynamics
Zhang et al. 3D inversion of large-scale frequency-domain airborne electromagnetic data using unstructured local mesh
CN117092702A (zh) 孔-隧激发极化探水结构的施工方法及反演探水方法
CN111158059A (zh) 基于三次b样条函数的重力反演方法
Mukanova et al. Modelling the influence of ground surface relief on electric sounding curves using the integral equations method

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
RJ01 Rejection of invention patent application after publication

Application publication date: 20210101

RJ01 Rejection of invention patent application after publication