CN113420487A - 一种三维重力梯度张量数值模拟方法、装置、设备和介质 - Google Patents

一种三维重力梯度张量数值模拟方法、装置、设备和介质 Download PDF

Info

Publication number
CN113420487A
CN113420487A CN202110978292.6A CN202110978292A CN113420487A CN 113420487 A CN113420487 A CN 113420487A CN 202110978292 A CN202110978292 A CN 202110978292A CN 113420487 A CN113420487 A CN 113420487A
Authority
CN
China
Prior art keywords
kernel function
function unit
dimensional
expression
calculated
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
CN202110978292.6A
Other languages
English (en)
Other versions
CN113420487B (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.)
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 CN202110978292.6A priority Critical patent/CN113420487B/zh
Publication of CN113420487A publication Critical patent/CN113420487A/zh
Application granted granted Critical
Publication of CN113420487B publication Critical patent/CN113420487B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • 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
    • G06F17/156—Correlation function computation including computation of convolution operations using a domain transform, e.g. Fourier transform, polynomial transform, number theoretic transform
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Mathematical Optimization (AREA)
  • Data Mining & Analysis (AREA)
  • General Engineering & Computer Science (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Computing Systems (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Algebra (AREA)
  • Geometry (AREA)
  • Databases & Information Systems (AREA)
  • Software Systems (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本申请涉及一种三维重力梯度张量数值模拟方法、装置、计算机设备和存储介质。所述方法包括：通过获取待计算区域在三维坐标系中的位置信息，将待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元，根据网格剖分和预先获取的核函数三维积分表达式得到核函数单元积分的解析表达式，根据核函数单元积分的解析表达式得到核函数单元系数矩阵，其中核函数单元系数矩阵为Toeplitz矩阵，通过三维离散快速傅里叶变换进行核函数单元系数矩阵和单元体密度值三维离散卷积的快速相乘算法，得到待计算区域的重力梯度张量异常值。本发明通过快速相乘方法，减小了系数矩阵存储和计算问题，大幅提高了正演的精度和效率。

Description

一种三维重力梯度张量数值模拟方法、装置、设备和介质
技术领域
本申请涉及计算机快速计算技术领域，特别是涉及一种三维重力梯度张量数值模拟方法、装置、计算机设备和存储介质。
背景技术
随着航空重力及测量仪器的不断进步， 重力梯度测量作为一种高精度、高分辨率进行地下目标精细勘探的重要手段，受到了国内众多学者的广泛关注。
重力梯度张量有5个独立分量，进行重力梯度测量对于地质解释而言，多分量组合反演解释能够得到更加合理的结果。当前，国内带重力梯度张量主要集中数据处理、数值模拟及反演解释等方面。其中，重力梯度张量数值模拟的研究相对较少，并且现有的方法计算较慢，无法实现大数据量快速地完成高精度反演解释工作。
针对上述问题，文献（Nagy,D., Papp, G., Benedek, J. The gravitationalpotential and its derivatives for the prism. Journal of Geodesy, 2000, 74(7-8):552-560）推导了棱柱体的引力位、重力异常及梯度张量的全空间解析表达式。文献(Pan, K. J, Zhang, Z. H, Hu, S. G, et al. Three-dimensional forward modellingof gravity field vector and its gradient tensor using the compact differencescheme. Geophysical Journal International, 2020, 224(2), 1272-1286) 采用高阶差分算子实现了基于有限差分法的重力场及梯度张量的正演计算，但随着网格单元数的增加，计算时间成倍的增长。（戴世坤, 陈轻蕊, 李昆,等. 重力异常场空间-波数混合域三维数值模拟. 地球物理学报, 2020, 063(005):2107-2119.）采用规则网格剖分，借助二维傅里叶变换实现了三维重力异常及梯度张量的快速数值模拟，该方法采用高斯傅里叶变换克服常规快速傅里叶变换的边界截断的影响，但计算异常体内部的梯度张量时误差相对较大，并且随着高斯点个数的增加计算量也会成倍增加。
已有的重力梯度张量三维数值模拟方法难以满足大数据量、快速三维反演解释的需求。因此，亟需研究一种高效、高精度的三维重力梯度张量数值模拟方法，以解决上述问题。
发明内容
基于此，有必要针对上述技术问题，提供一种能够提高三维重力梯度张量数值模拟效率和精度的三维重力梯度张量数值模拟方法、装置、计算机设备和存储介质。
一种三维重力梯度张量数值模拟方法，所述方法包括：
获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；
根据待计算区域中三维地质体的形状和密度设定所述网格单元的密度，得到单元体密度值；
通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值。
在其中一个实施例中，还包括：构建三维坐标系，其中x轴指向正东，y轴指向正北，z轴垂直向下；
获取待计算区域在所述三维坐标系中的位置信息；
将所述待计算区域沿x、y、z方向分别等间隔剖分成
Figure 489898DEST_PATH_IMAGE001
个小长方体，网格单元间隔分别为Δx、Δy、Δz，形成多个长方体的网格单元。
在其中一个实施例中，还包括：以所述网格单元的中心点位置为观测点；
获取重力梯度垂向分量核函数三维积分表达式为：
Figure 883970DEST_PATH_IMAGE002
其中，G表示引力常数，
Figure 577120DEST_PATH_IMAGE004
为观测点坐标，
Figure 107458DEST_PATH_IMAGE005
表示源点坐标；
对所述重力梯度垂向分量核函数三维积分表达式进行离散，得到重力梯度垂向分量核函数单元离散表达式为：
Figure 657388DEST_PATH_IMAGE006
其中，
Figure 273177DEST_PATH_IMAGE007
表示编号为(i, j, k)的网格单元的中心坐标，
Figure 566756DEST_PATH_IMAGE008
表示编号为(l, m, n)的网格单元的中心坐标,其中，
Figure 217180DEST_PATH_IMAGE009
；
Figure 108650DEST_PATH_IMAGE010
；
Figure 477315DEST_PATH_IMAGE011
；编号为(l, m, n)的网格单元三个方向积分区间分别为
Figure 309004DEST_PATH_IMAGE012
至
Figure 79514DEST_PATH_IMAGE014
、
Figure 908930DEST_PATH_IMAGE016
至
Figure 764890DEST_PATH_IMAGE017
、
Figure 337954DEST_PATH_IMAGE018
至
Figure 697391DEST_PATH_IMAGE019
。
在其中一个实施例中，还包括：对所述重力梯度垂向分量核函数单元离散表达式进行积分，得到重力梯度垂向分量核函数单元积分的解析表达式为：
Figure 963288DEST_PATH_IMAGE020
其中，
Figure 837703DEST_PATH_IMAGE021
为整数，
Figure 712993DEST_PATH_IMAGE022
，
Figure 192516DEST_PATH_IMAGE023
，
Figure 426051DEST_PATH_IMAGE024
、
Figure 991024DEST_PATH_IMAGE025
、
Figure 905891DEST_PATH_IMAGE026
、
Figure 239920DEST_PATH_IMAGE028
、
Figure 644357DEST_PATH_IMAGE029
，
Figure 431047DEST_PATH_IMAGE030
。
在其中一个实施例中，还包括：根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；所述核函数单元系数矩阵大小为
Figure 149604DEST_PATH_IMAGE031
。
在其中一个实施例中，还包括：通过三维离散快速傅里叶变换将所述核函数单元系数矩阵和所述单元体密度值的三维离散卷积计算转换为快速相乘算法，得到所述待计算区域的重力梯度张量为：
Figure 603719DEST_PATH_IMAGE032
其中，F和
Figure 880855DEST_PATH_IMAGE033
分别表示三维离散傅里叶正反变换算子；
Figure 686000DEST_PATH_IMAGE034
为源项密度扩展矩阵，
Figure 208248DEST_PATH_IMAGE036
为单元体密度矩阵，
Figure 251290DEST_PATH_IMAGE037
为扩展矩阵，
Figure 997529DEST_PATH_IMAGE038
表示提取矩阵的前
Figure 493233DEST_PATH_IMAGE039
个元素；
根据计算结果得到所述待计算区域的重力梯度张量异常值。
在其中一个实施例中，还包括：不同网格单元的单元体密度值不同。
一种三维重力梯度张量数值模拟装置，所述装置包括：
网格剖分模块，用于获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
核函数单元系数矩阵确定模块，用于根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；
单元体密度值确定模块，用于根据待计算区域中三维地质体的形状和密度设定所述网格单元的密度，得到单元体密度值；
重力梯度张量异常值计算模块，用于通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值。
一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，所述处理器执行所述计算机程序时实现以下步骤：
获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；
根据待计算区域中三维地质体的形状和密度设定所述网格单元的密度，得到单元体密度值；
通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值。
一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执行时实现以下步骤：
获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；
根据待计算区域中三维地质体的形状和密度设定所述网格单元的密度，得到单元体密度值；
通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值。
上述三维重力梯度张量数值模拟方法、装置、计算机设备和存储介质，通过获取待计算区域在三维坐标系中的位置信息，将待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据核函数单元积分的解析表达式得到核函数单元系数矩阵，其中核函数单元系数矩阵为Toeplitz矩阵；根据待计算区域中三维地质体的形状和密度设定网格单元的密度，得到单元体密度值；通过三维离散快速傅里叶变换进行核函数单元系数矩阵和单元体密度值三维离散卷积的快速相乘算法，得到待计算区域的重力梯度张量异常值。本发明每个网格单元的核函数积分都可推导出解析解，提高了计算精度，充分利用Toeplitz型矩阵的特点，大大降低核函数矩阵的存储空间和计算时间；同时，借助三维离散快速傅里叶变换实现核函数单元积分系数矩阵与密度快速相乘，一次可计算整个三维网格单元上（包括异常体内部）的重力梯度张量异常值，提高了正演计算的效率。
附图说明
图1为一个实施例中三维重力梯度张量数值模拟方法的流程场景图；
图2为一个实施例中计算区域网格剖分示意图；
图3为一个具体实施例中使用本发明方法计算的重力梯度垂向分量三维数值模拟结果图；
图4为一个具体实施例中解析解计算的重力梯度垂向分量三维结果图；
图5为一个具体实施例中重力梯度垂向分量使用本发明方法的数值模拟结果与解析解的相对误差三维结果图；
图6为一个具体实施例中使用本发明方法计算的地表重力梯度垂向分量图；
图7为一个具体实施例中解析解计算的地表重力梯度垂向分量图
图8为一个具体实施例中地表重力梯度垂向分量使用本发明方法的数值模拟结果与解析解的相对误差三维结果图；
图9为一个实施例中三维重力梯度张量数值模拟装置的结构框图；
图10为一个实施例中计算机设备的内部结构图。
具体实施方式
为了使本申请的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本申请进行进一步详细说明。应当理解，此处描述的具体实施例仅仅用以解释本申请，并不用于限定本申请。
本申请提供的三维重力梯度张量数值模拟方法，可以应用于如下应用环境中。其中，终端执行一种三维重力梯度张量数值模拟方法，通过获取待计算区域在三维坐标系中的位置信息，将待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据核函数单元积分的解析表达式得到核函数单元系数矩阵，其中核函数单元系数矩阵为Toeplitz矩阵；根据待计算区域中三维地质体的形状和密度设定网格单元的密度，得到单元体密度值；通过三维离散快速傅里叶变换进行核函数单元系数矩阵和单元体密度值三维离散卷积的快速相乘算法，得到待计算区域的重力梯度张量异常值。其中，终端可以但不限于是各种个人计算机、笔记本电脑、平板电脑和便携式设备。
在一个实施例中，如图1所示，提供了一种三维重力梯度张量数值模拟方法，包括以下步骤：
步骤102，获取待计算区域在三维坐标系中的位置信息，将待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元。
如图2所示，将该计算区域沿着x、y、z方向分别等间隔剖分成
Figure 819172DEST_PATH_IMAGE040
个小长方体，网格单元间隔分别为Δx、Δy、Δz。观测点的位置和每个单元的中心点位置重合，观测点的总数为
Figure 982300DEST_PATH_IMAGE041
个。
步骤104，根据网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据核函数单元积分的解析表达式得到核函数单元系数矩阵。
地下三维地质体在观测点上产生的重力梯度张量可以表示为核函数与密度的卷积，例如重力梯度垂向分量为：
Figure 899440DEST_PATH_IMAGE042
式中，
Figure 413598DEST_PATH_IMAGE043
表示重力异常核函数单元积分系数矩阵，
Figure 776184DEST_PATH_IMAGE045
表示三维卷积算子，
Figure 793819DEST_PATH_IMAGE047
表示每个剖分单元的密度值，
Figure 85123DEST_PATH_IMAGE048
表示重力梯度垂向分量。
核函数单元系数矩阵为Toeplitz矩阵，Toeplitz矩阵只需计算和存储第一行和第一列
Figure 820998DEST_PATH_IMAGE049
，进一步将t写成
Figure 285477DEST_PATH_IMAGE050
，则
Figure 892039DEST_PATH_IMAGE051
。对称的Toeplitz型矩阵，可以大大降低核函数矩阵的存储空间和计算时间。
步骤106，根据待计算区域中三维地质体的形状和密度设定网格单元的密度，得到单元体密度值。
根据三维地质体的形状和密度分布设定网格剖分的每个单元的密度值；每个单元体的密度值为常值，不同网格单元体的密度值不同，以此描述任意密度分布的复杂地质体模型。
步骤108，通过三维离散快速傅里叶变换进行核函数单元系数矩阵和单元体密度值三维离散卷积的快速相乘算法，得到待计算区域的重力梯度张量异常值。
当通过一维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值一维离散卷积的快速相乘算法时，一维离散卷积可以表示为：
Figure 619824DEST_PATH_IMAGE052
根据所述重力梯度垂向分量核函数单元积分的解析表达式得到核函数单元系数矩阵为：
Figure 639732DEST_PATH_IMAGE054
进一步表示为：
Figure 580006DEST_PATH_IMAGE055
其中，核函数单元系数矩阵
Figure 572233DEST_PATH_IMAGE056
为Toeplitz矩阵，令
Figure 267657DEST_PATH_IMAGE057
，则
Figure 211080DEST_PATH_IMAGE058
，进一步地，令
Figure 955045DEST_PATH_IMAGE059
，则
Figure 67357DEST_PATH_IMAGE060
，
Figure 871365DEST_PATH_IMAGE061
表示扩展矩阵；
所述核函数单元系数矩阵和所述单元体密度值的一维离散卷积可以表示为Toeplitz矩阵
Figure 69128DEST_PATH_IMAGE061
与向量的快速相乘，为：
Figure 413522DEST_PATH_IMAGE062
其中，
Figure 114762DEST_PATH_IMAGE063
和
Figure 89671DEST_PATH_IMAGE064
分别表示一维离散傅里叶正反变换算子；一维源向量
Figure 571468DEST_PATH_IMAGE065
，
Figure 657236DEST_PATH_IMAGE066
表示提取矩阵的前
Figure 275299DEST_PATH_IMAGE067
个元素。
三维离散卷积核函数单元积分系数矩阵与密度快速相乘与一维离散卷积类似。通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法可以一次计算获得整个三维剖分网格单元（包括异常体内部）上的重力梯度张量异常值。
通过系数矩阵的对称性和核函数单元积分系数矩阵与密度的快速相乘方法，减小了系数矩阵存储和计算问题，大幅提高了正演的精度和效率。
上述三维重力梯度张量数值模拟方法中，通过获取待计算区域在三维坐标系中的位置信息，将待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据核函数单元积分的解析表达式得到核函数单元系数矩阵，其中核函数单元系数矩阵为Toeplitz矩阵；根据待计算区域中三维地质体的形状和密度设定网格单元的密度，得到单元体密度值；通过三维离散快速傅里叶变换进行核函数单元系数矩阵和单元体密度值三维离散卷积的快速相乘算法，得到待计算区域的重力梯度张量异常值。本发明每个网格单元的核函数积分都可推导出解析解，提高了计算精度，充分利用Toeplitz型矩阵的特点，大大降低核函数矩阵的存储空间和计算时间；同时，借助三维离散快速傅里叶变换实现核函数单元积分系数矩阵与密度快速相乘，一次可计算整个三维网格单元上（包括异常体内部）的重力梯度张量异常值，提高了正演计算的效率。
在其中一个实施例中，还包括：构建三维坐标系，其中x轴指向正东，y轴指向正北，z轴垂直向下；获取待计算区域在三维坐标系中的位置信息；将待计算区域沿x、y、z方向分别等间隔剖分成
Figure 185224DEST_PATH_IMAGE068
个小长方体，网格单元间隔分别为Δx、Δy、Δz，形成多个长方体的网格单元。
在其中一个实施例中，还包括：以网格单元的中心点位置为观测点；
获取重力梯度垂向分量核函数三维积分表达式为：
Figure 826421DEST_PATH_IMAGE069
其中，G表示引力常数，
Figure 512617DEST_PATH_IMAGE070
为观测点坐标，
Figure 454028DEST_PATH_IMAGE071
表示源点坐标；
对重力梯度垂向分量核函数三维积分表达式进行离散，得到重力梯度垂向分量核函数单元离散表达式为：
Figure 505161DEST_PATH_IMAGE072
其中，
Figure 899233DEST_PATH_IMAGE073
表示编号为(i, j, k)的网格单元的中心坐标，
Figure 389120DEST_PATH_IMAGE074
表示编号为(l, m, n)的网格单元的中心坐标,其中，
Figure 919459DEST_PATH_IMAGE075
；
Figure 672651DEST_PATH_IMAGE076
；
Figure 350757DEST_PATH_IMAGE077
；编号为(l, m, n)的网格单元三个方向积分区间分别为
Figure 316439DEST_PATH_IMAGE079
至
Figure 465398DEST_PATH_IMAGE080
、
Figure 186230DEST_PATH_IMAGE081
至
Figure 289315DEST_PATH_IMAGE082
、
Figure 324267DEST_PATH_IMAGE084
至
Figure 625935DEST_PATH_IMAGE085
。
在其中一个实施例中，还包括：对重力梯度垂向分量核函数单元离散表达式进行积分，得到重力梯度垂向分量核函数单元积分的解析表达式为：
Figure 720930DEST_PATH_IMAGE086
其中，
Figure 311312DEST_PATH_IMAGE087
为整数，
Figure 884375DEST_PATH_IMAGE088
，
Figure 306129DEST_PATH_IMAGE089
，
Figure 572026DEST_PATH_IMAGE090
、
Figure 384124DEST_PATH_IMAGE091
、
Figure 823195DEST_PATH_IMAGE092
、
Figure 801253DEST_PATH_IMAGE093
、
Figure 972472DEST_PATH_IMAGE094
，
Figure 334183DEST_PATH_IMAGE095
。
每个单元的核函数积分都可以推导出解析解，进一步提高计算精度。
在其中一个实施例中，还包括：通过三维离散快速傅里叶变换将核函数单元系数矩阵和单元体密度值的三维离散卷积计算转换为快速相乘算法，得到待计算区域的重力梯度张量为：
Figure 983470DEST_PATH_IMAGE096
其中，F和
Figure 583079DEST_PATH_IMAGE097
分别表示三维离散傅里叶正反变换算子；
Figure 987515DEST_PATH_IMAGE098
为源项密度扩展矩阵，
Figure 39785DEST_PATH_IMAGE099
为单元体密度矩阵，
Figure 492763DEST_PATH_IMAGE100
为扩展矩阵，
Figure 946878DEST_PATH_IMAGE101
表示提取矩阵的前
Figure 522216DEST_PATH_IMAGE102
个元素；
根据计算结果得到待计算区域的重力梯度张量异常值。
以上实施例以计算重力梯度垂向分量为例进行分析，别的分量也可以按照该算法进行计算，只需改变梯度分量满足的三维积分表达式。
在其中一个实施例中，还包括：不同网格单元的单元体密度值不同。
应该理解的是，虽然图1的流程图中的各个步骤按照箭头的指示依次显示，但是这些步骤并不是必然按照箭头指示的顺序依次执行。除非本文中有明确的说明，这些步骤的执行并没有严格的顺序限制，这些步骤可以以其它的顺序执行。而且，图1中的至少一部分步骤可以包括多个子步骤或者多个阶段，这些子步骤或者阶段并不必然是在同一时刻执行完成，而是可以在不同的时刻执行，这些子步骤或者阶段的执行顺序也不必然是依次进行，而是可以与其它步骤或者其它步骤的子步骤或者阶段的至少一部分轮流或者交替地执行。
在一个具体实施例中，计算区域内有一个棱柱型异常体，计算区域范围为：x和y方向均从-500m到500m，z方向从0m到1000m(z轴垂直向下为正)，网格剖分间隔均为10m，整个计算区域剖分为100×100×100个单元，异常体分布范围为：x和y方向均从-50m到50m，z方向从150m到250m，剩余密度为2000kg/m3，计算整个三维网格上的重力异常。
利用Fortran语言编程实现，运行程序所用的个人电脑配置为：CPU-Inter Corei7-8700，主频为3.2GHz。计算整个网格单元上的值，即1×106个观测点的值需要0.558秒，可以看出该方法计算效率很高。图3、图4分别为本发明方法和解析解计算的重力梯度垂向分量三维结果图，从两幅图可以看出解析解与本发明方法计算结果形态吻合很好，通过对三维立体图分析可以看出重力梯度垂向分量能够很好反演异常的边界；图5为重力梯度垂向分量的解析解和本发明方法计算的相对误差，可以看出整个三维网格单元上，即异常内部和外部的相对误差都小于7.1×10-9，进一步可以看出该方法对于异常体内部和外部都可以达到很高的计算精度；图6、图7分别为图3和图4中抽取地表的本发明方法和解析解计算的重力梯度垂向分量，可以看出两幅图形态一致，图8为地面的重力梯度垂向分量相对误差图，可以看出整个平面相对误差都小于4.0×10-10；表1给出了重力梯度垂向分量地面的相对误差进行统计，从表1可以看出本方法计算精度较高。
表1本发明方法计算值与重力梯度垂向分量理论值的误差统计表（单位：mGal）
Figure 61782DEST_PATH_IMAGE103
在一个实施例中，如图9所示，提供了一种三维重力梯度张量数值模拟装置，包括：网格剖分模块902、核函数单元系数矩阵确定模块904、单元体密度值确定模块906和重力梯度张量异常值计算模块908，其中：
网格剖分模块902，用于获取待计算区域在三维坐标系中的位置信息，将待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
核函数单元系数矩阵确定模块904，用于根据网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据核函数单元积分的解析表达式得到核函数单元系数矩阵；核函数单元系数矩阵为Toeplitz矩阵；
单元体密度值确定模块906，用于根据待计算区域中三维地质体的形状和密度设定网格单元的密度，得到单元体密度值；
重力梯度张量异常值计算模块908，用于通过三维离散快速傅里叶变换进行核函数单元系数矩阵和单元体密度值三维离散卷积的快速相乘算法，得到待计算区域的重力梯度张量异常值。
网格剖分模块902还用于构建三维坐标系，其中x轴指向正东，y轴指向正北，z轴垂直向下；
获取待计算区域在三维坐标系中的位置信息；
将待计算区域沿x、y、z方向分别等间隔剖分成
Figure 816986DEST_PATH_IMAGE104
个小长方体，网格单元间隔分别为Δx、Δy、Δz，形成多个长方体的网格单元。
核函数单元系数矩阵确定模块904还用于以网格单元的中心点位置为观测点；
获取重力梯度垂向分量核函数三维积分表达式为：
Figure 125607DEST_PATH_IMAGE105
其中，G表示引力常数，
Figure 871847DEST_PATH_IMAGE107
为观测点坐标，
Figure 367550DEST_PATH_IMAGE108
表示源点坐标；
对重力梯度垂向分量核函数三维积分表达式进行离散，得到重力梯度垂向分量核函数单元离散表达式为：
Figure 693489DEST_PATH_IMAGE109
其中，
Figure 653355DEST_PATH_IMAGE110
表示编号为(i, j, k)的网格单元的中心坐标，
Figure 508178DEST_PATH_IMAGE112
表示编号为(l, m, n)的网格单元的中心坐标,其中，
Figure 22336DEST_PATH_IMAGE113
；
Figure 886387DEST_PATH_IMAGE114
；
Figure 700759DEST_PATH_IMAGE115
；编号为(l, m, n)的网格单元三个方向积分区间分别为
Figure 992063DEST_PATH_IMAGE117
至
Figure 226473DEST_PATH_IMAGE118
、
Figure 425374DEST_PATH_IMAGE119
至
Figure 563094DEST_PATH_IMAGE120
、
Figure 25299DEST_PATH_IMAGE121
至
Figure 779629DEST_PATH_IMAGE122
。
核函数单元系数矩阵确定模块904还用于对重力梯度垂向分量核函数单元离散表达式进行积分，得到重力梯度垂向分量核函数单元积分的解析表达式为：
Figure 251061DEST_PATH_IMAGE123
其中，
Figure 243288DEST_PATH_IMAGE124
为整数，
Figure 610815DEST_PATH_IMAGE125
，
Figure 321282DEST_PATH_IMAGE126
，
Figure 127564DEST_PATH_IMAGE127
、
Figure 708718DEST_PATH_IMAGE128
、
Figure 575043DEST_PATH_IMAGE129
、
Figure 5762DEST_PATH_IMAGE130
、
Figure 553418DEST_PATH_IMAGE131
，
Figure DEST_PATH_IMAGE132
。
重力梯度张量异常值计算模块908还用于通过三维离散快速傅里叶变换将核函数单元系数矩阵和单元体密度值的三维离散卷积计算转换为快速相乘算法，得到待计算区域的重力梯度张量为：
Figure 785817DEST_PATH_IMAGE133
其中，F和
Figure DEST_PATH_IMAGE134
分别表示三维离散傅里叶正反变换算子；
Figure 963988DEST_PATH_IMAGE135
为源项密度扩展矩阵，
Figure 445785DEST_PATH_IMAGE137
为单元体密度矩阵，
Figure DEST_PATH_IMAGE138
为扩展矩阵，
Figure 394DEST_PATH_IMAGE139
表示提取矩阵的前
Figure DEST_PATH_IMAGE140
个元素；
根据计算结果得到所述待计算区域的重力梯度张量异常值。
关于三维重力梯度张量数值模拟装置的具体限定可以参见上文中对于三维重力梯度张量数值模拟方法的限定，在此不再赘述。上述三维重力梯度张量数值模拟装置中的各个模块可全部或部分通过软件、硬件及其组合来实现。上述各模块可以硬件形式内嵌于或独立于计算机设备中的处理器中，也可以以软件形式存储于计算机设备中的存储器中，以便于处理器调用执行以上各个模块对应的操作。
在一个实施例中，提供了一种计算机设备，该计算机设备可以是终端，其内部结构图可以如图10所示。该计算机设备包括通过系统总线连接的处理器、存储器、网络接口、显示屏和输入装置。其中，该计算机设备的处理器用于提供计算和控制能力。该计算机设备的存储器包括非易失性存储介质、内存储器。该非易失性存储介质存储有操作系统和计算机程序。该内存储器为非易失性存储介质中的操作系统和计算机程序的运行提供环境。该计算机设备的网络接口用于与外部的终端通过网络连接通信。该计算机程序被处理器执行时以实现一种三维重力梯度张量数值模拟方法。该计算机设备的显示屏可以是液晶显示屏或者电子墨水显示屏，该计算机设备的输入装置可以是显示屏上覆盖的触摸层，也可以是计算机设备外壳上设置的按键、轨迹球或触控板，还可以是外接的键盘、触控板或鼠标等。
本领域技术人员可以理解，图10中示出的结构，仅仅是与本申请方案相关的部分结构的框图，并不构成对本申请方案所应用于其上的计算机设备的限定，具体的计算机设备可以包括比图中所示更多或更少的部件，或者组合某些部件，或者具有不同的部件布置。
在一个实施例中，提供了一种计算机设备，包括存储器和处理器，该存储器存储有计算机程序，该处理器执行计算机程序时实现上述方法实施例中的步骤。
在一个实施例中，提供了一种计算机可读存储介质，其上存储有计算机程序，计算机程序被处理器执行时实现上述方法实施例中的步骤。
本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程，是可以通过计算机程序来指令相关的硬件来完成，所述的计算机程序可存储于一非易失性计算机可读取存储介质中，该计算机程序在执行时，可包括如上述各方法的实施例的流程。其中，本申请所提供的各实施例中所使用的对存储器、存储、数据库或其它介质的任何引用，均可包括非易失性和/或易失性存储器。非易失性存储器可包括只读存储器（ROM）、可编程ROM（PROM）、电可编程ROM（EPROM）、电可擦除可编程ROM（EEPROM）或闪存。易失性存储器可包括随机存取存储器（RAM）或者外部高速缓冲存储器。作为说明而非局限，RAM以多种形式可得，诸如静态RAM（SRAM）、动态RAM（DRAM）、同步DRAM（SDRAM）、双数据率SDRAM（DDRSDRAM）、增强型SDRAM（ESDRAM）、同步链路（Synchlink） DRAM（SLDRAM）、存储器总线（Rambus）直接RAM（RDRAM）、直接存储器总线动态RAM（DRDRAM）、以及存储器总线动态RAM（RDRAM）等。
以上实施例的各技术特征可以进行任意的组合，为使描述简洁，未对上述实施例中的各个技术特征所有可能的组合都进行描述，然而，只要这些技术特征的组合不存在矛盾，都应当认为是本说明书记载的范围。
以上所述实施例仅表达了本申请的几种实施方式，其描述较为具体和详细，但并不能因此而理解为对发明专利范围的限制。应当指出的是，对于本领域的普通技术人员来说，在不脱离本申请构思的前提下，还可以做出若干变形和改进，这些都属于本申请的保护范围。因此，本申请专利的保护范围应以所附权利要求为准。

Claims (10)

1.一种三维重力梯度张量数值模拟方法，其特征在于，所述方法包括：
获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；
根据待计算区域中三维地质体的形状和密度设定所述网格单元的密度，得到单元体密度值；
通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值。
2.根据权利要求1所述的方法，其特征在于，获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元，包括：
构建三维坐标系，其中x轴指向正东，y轴指向正北，z轴垂直向下；
获取待计算区域在所述三维坐标系中的位置信息；
将所述待计算区域沿x、y、z方向分别等间隔剖分成
Figure 119534DEST_PATH_IMAGE001
个小长方体，网格单元间隔分别为Δx、Δy、Δz，形成多个长方体的网格单元。
3.根据权利要求2所述的方法，其特征在于，根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，包括：
以所述网格单元的中心点位置为观测点；
获取重力梯度垂向分量核函数三维积分表达式为：
Figure 505516DEST_PATH_IMAGE002
其中，G表示引力常数，
Figure 488515DEST_PATH_IMAGE004
为观测点坐标，
Figure 149304DEST_PATH_IMAGE005
表示源点坐标；
对所述重力梯度垂向分量核函数三维积分表达式进行离散，得到重力梯度垂向分量核函数单元离散表达式为：
Figure 432517DEST_PATH_IMAGE006
其中，
Figure 956777DEST_PATH_IMAGE007
表示编号为(i, j, k)的网格单元的中心坐标，
Figure 489390DEST_PATH_IMAGE009
表示编号为(l, m, n)的网格单元的中心坐标,其中，
Figure 157132DEST_PATH_IMAGE010
；
Figure 29273DEST_PATH_IMAGE011
；
Figure 491478DEST_PATH_IMAGE012
；编号为(l, m, n)的网格单元三个方向积分区间分别为
Figure 511387DEST_PATH_IMAGE013
至
Figure 451661DEST_PATH_IMAGE015
、
Figure 443888DEST_PATH_IMAGE017
至
Figure 139311DEST_PATH_IMAGE018
、
Figure 584199DEST_PATH_IMAGE019
至
Figure 328164DEST_PATH_IMAGE020
。
4.根据权利要求3所述的方法，其特征在于，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，包括：
对所述重力梯度垂向分量核函数单元离散表达式进行积分，得到重力梯度垂向分量核函数单元积分的解析表达式为：
Figure 237214DEST_PATH_IMAGE021
其中，
Figure 274178DEST_PATH_IMAGE023
为整数，
Figure 471941DEST_PATH_IMAGE024
，
Figure 816335DEST_PATH_IMAGE025
，
Figure 517575DEST_PATH_IMAGE026
、
Figure 758063DEST_PATH_IMAGE027
、
Figure 177543DEST_PATH_IMAGE028
、
Figure 325628DEST_PATH_IMAGE029
、
Figure 146953DEST_PATH_IMAGE030
，
Figure 292764DEST_PATH_IMAGE031
。
5.根据权利要求4所述的方法，其特征在于，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵，包括：
根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；所述核函数单元系数矩阵大小为
Figure 199540DEST_PATH_IMAGE032
。
6.根据权利要求5所述的方法，其特征在于，通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值，包括：
通过三维离散快速傅里叶变换将所述核函数单元系数矩阵和所述单元体密度值的三维离散卷积计算转换为快速相乘算法，得到所述待计算区域的重力梯度张量为：
Figure 885736DEST_PATH_IMAGE033
其中，F和
Figure 60103DEST_PATH_IMAGE034
分别表示三维离散傅里叶正反变换算子；
Figure DEST_PATH_IMAGE035
为源项密度扩展矩阵，
Figure DEST_PATH_IMAGE037
为单元体密度矩阵，
Figure DEST_PATH_IMAGE038
为扩展矩阵，
Figure DEST_PATH_IMAGE039
表示提取矩阵的前
Figure DEST_PATH_IMAGE040
个元素；
根据计算结果得到所述待计算区域的重力梯度张量异常值。
7.根据权利要求1至6任意一项所述的方法，其特征在于，不同网格单元的单元体密度值不同。
8.一种三维重力梯度张量数值模拟装置，其特征在于，所述装置包括：
网格剖分模块，用于获取待计算区域在三维坐标系中的位置信息，将所述待计算区域沿x、y、z方向分别进行等间隔剖分，形成多个长方体的网格单元；
核函数单元系数矩阵确定模块，用于根据所述网格单元，对预先获取的核函数三维积分表达式进行离散，得到核函数单元离散表达式，对所述核函数单元离散表达式进行积分，得到核函数单元积分的解析表达式，根据所述核函数单元积分的解析表达式得到核函数单元系数矩阵；所述核函数单元系数矩阵为Toeplitz矩阵；
单元体密度值确定模块，用于根据待计算区域中三维地质体的形状和密度设定所述网格单元的密度，得到单元体密度值；
重力梯度张量异常值计算模块，用于通过三维离散快速傅里叶变换进行所述核函数单元系数矩阵和所述单元体密度值三维离散卷积的快速相乘算法，得到所述待计算区域的重力梯度张量异常值。
9.一种计算机设备，包括存储器和处理器，所述存储器存储有计算机程序，其特征在于，所述处理器执行所述计算机程序时实现权利要求1至7中任一项所述方法的步骤。
10.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述计算机程序被处理器执行时实现权利要求1至7中任一项所述的方法的步骤。
CN202110978292.6A 2021-08-25 2021-08-25 一种三维重力梯度张量数值模拟方法、装置、设备和介质 Active CN113420487B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110978292.6A CN113420487B (zh) 2021-08-25 2021-08-25 一种三维重力梯度张量数值模拟方法、装置、设备和介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110978292.6A CN113420487B (zh) 2021-08-25 2021-08-25 一种三维重力梯度张量数值模拟方法、装置、设备和介质

Publications (2)

Publication Number Publication Date
CN113420487A true CN113420487A (zh) 2021-09-21
CN113420487B CN113420487B (zh) 2021-10-29

Family

ID=77719386

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110978292.6A Active CN113420487B (zh) 2021-08-25 2021-08-25 一种三维重力梯度张量数值模拟方法、装置、设备和介质

Country Status (1)

Country Link
CN (1) CN113420487B (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113806686A (zh) * 2021-11-19 2021-12-17 中南大学 大规模复杂地质体重力梯度快速计算方法、装置和设备
CN116609844A (zh) * 2023-03-13 2023-08-18 中国自然资源航空物探遥感中心 一种用于提取地球物理三维断裂信息的方法和计算设备
CN116911146A (zh) * 2023-09-14 2023-10-20 中南大学 三维重力场全息数值模拟及cpu-gpu加速方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20110015909A1 (en) * 2009-07-16 2011-01-20 Gang Zhao Reservoir modeling method
CN103163562A (zh) * 2013-02-01 2013-06-19 中国科学院测量与地球物理研究所 基于滤波原理的卫星重力梯度反演方法
CN107577641A (zh) * 2017-08-21 2018-01-12 吉林大学 一种基于gpu并行的重力梯度张量数据快速密度反演方法
CN109375280A (zh) * 2018-12-10 2019-02-22 中南大学 一种球坐标系下重力场快速高精度正演方法
US20200066035A1 (en) * 2018-08-27 2020-02-27 Exxonmobil Research And Engineering Company Inversion of large, nearly-homogeneous geobodies via ultra low-dimensional shape representation using orthogonal basis functions
CN111400654A (zh) * 2020-03-13 2020-07-10 中南大学 一种基于Toplitze核矩阵的重力场快速正演方法及反演方法
CN112800657A (zh) * 2021-04-15 2021-05-14 中南大学 基于复杂地形的重力场数值模拟方法、装置和计算机设备

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20110015909A1 (en) * 2009-07-16 2011-01-20 Gang Zhao Reservoir modeling method
CN103163562A (zh) * 2013-02-01 2013-06-19 中国科学院测量与地球物理研究所 基于滤波原理的卫星重力梯度反演方法
CN107577641A (zh) * 2017-08-21 2018-01-12 吉林大学 一种基于gpu并行的重力梯度张量数据快速密度反演方法
US20200066035A1 (en) * 2018-08-27 2020-02-27 Exxonmobil Research And Engineering Company Inversion of large, nearly-homogeneous geobodies via ultra low-dimensional shape representation using orthogonal basis functions
CN109375280A (zh) * 2018-12-10 2019-02-22 中南大学 一种球坐标系下重力场快速高精度正演方法
CN111400654A (zh) * 2020-03-13 2020-07-10 中南大学 一种基于Toplitze核矩阵的重力场快速正演方法及反演方法
CN112800657A (zh) * 2021-04-15 2021-05-14 中南大学 基于复杂地形的重力场数值模拟方法、装置和计算机设备

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
XULONG WANG: "Efficient 2D Modeling of Gravity Anomaly in Space-wavenumber Mixed Domain", 《9TH INTERNATIONAL CONFERENCE ON ENVIRONMENTAL AND ENGINEERING GEOPHYSICS》 *
张琦: "重力梯度数据的边界识别与深度成像方法研究", 《中国优秀硕士学位论文全文数据库 基础科学辑》 *
王旭龙等: "空间—波数混合域二度体重力异常正演", 《石油地球物理勘探》 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113806686A (zh) * 2021-11-19 2021-12-17 中南大学 大规模复杂地质体重力梯度快速计算方法、装置和设备
CN113806686B (zh) * 2021-11-19 2022-02-08 中南大学 大规模复杂地质体重力梯度快速计算方法、装置和设备
CN116609844A (zh) * 2023-03-13 2023-08-18 中国自然资源航空物探遥感中心 一种用于提取地球物理三维断裂信息的方法和计算设备
CN116609844B (zh) * 2023-03-13 2024-03-12 中国自然资源航空物探遥感中心 一种用于提取地球物理三维断裂信息的方法和计算设备
CN116911146A (zh) * 2023-09-14 2023-10-20 中南大学 三维重力场全息数值模拟及cpu-gpu加速方法
CN116911146B (zh) * 2023-09-14 2024-01-19 中南大学 三维重力场全息数值模拟及cpu-gpu加速方法

Also Published As

Publication number Publication date
CN113420487B (zh) 2021-10-29

Similar Documents

Publication Publication Date Title
CN113420487B (zh) 一种三维重力梯度张量数值模拟方法、装置、设备和介质
CN112800657B (zh) 基于复杂地形的重力场数值模拟方法、装置和计算机设备
Spitzer A 3-D finite-difference algorithm for DC resistivity modelling using conjugate gradient methods
CN112287534B (zh) 基于nufft的二维磁异常快速正演模拟方法和装置
Li et al. A space–time fractional phase-field model with tunable sharpness and decay behavior and its efficient numerical simulation
Jablonowski et al. Block-structured adaptive grids on the sphere: Advection experiments
Shankar et al. Mesh-free semi-Lagrangian methods for transport on a sphere using radial basis functions
CN111967169B (zh) 二度体重力异常积分解数值模拟方法和装置
Zhang et al. Meshless regular hybrid boundary node method
Kruglyakov et al. High-performance parallel solver for integral equations of electromagnetics based on Galerkin method
Zhang et al. A hybrid iterative numerical transferable solver (HINTS) for PDEs based on deep operator network and relaxation methods
Barral et al. Anisotropic mesh adaptation in Firedrake with PETSc DMPlex
Hu et al. A meshless collocation method for solving the inverse Cauchy problem associated with the variable-order fractional heat conduction model under functionally graded materials
CN113673163B (zh) 一种三维磁异常数快速正演方法、装置和计算机设备
Chen et al. Numerical approximation of the inviscid 3D primitive equations in a limited domain
CN113076678B (zh) 一种频率域二度体重力异常快速数值模拟方法和装置
CN113221409B (zh) 一种有限元和边界元耦合的声波二维数值模拟方法和装置
Siebenborn et al. A curved-element unstructured discontinuous Galerkin method on GPUs for the Euler equations
Yang et al. Parallel domain decomposition methods with mixed order discretization for fully implicit solution of tracer transport problems on the cubed-sphere
CN113642189B (zh) 一种基于积分解的重力梯度张量快速数值模拟方法和装置
CN113806686B (zh) 大规模复杂地质体重力梯度快速计算方法、装置和设备
Tibaut et al. Fast Boundary-Domain Integral Method with the H2-matrix formulation for large scale numerical investigations
CN113656976A (zh) 一种二维磁梯度张量快速数值模拟方法、装置和设备
Šprlák et al. Vertical and horizontal spheroidal boundary-value problems
Ostermann Three-dimensional modeling of heat transport in deep hydrothermal reservoirs

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