三维复杂地质结构模型的InterfaceGrid表达方法

InterfaceGrid: Gridding representation of 3D geological models for complex geological structures

摘要

三维地质结构模型是地质空间内地质构造和地质体对象边界的数字化表征模型。随着地质数据矢量栅格一体化集成和地质数据空间查询分析需求增多,构建统一表达地质结构的一体化空间数据模型成为三维地质信息科学领域需要解决的基础性问题之一。针对规则网格、PillarGrid和SBRT模型等难以实现复杂地质结构的精细表达问题,本研究提出InterfaceGrid数据模型,充分考虑地质结构和属性的地下分布具有极强的非均匀性、不连续性和空间多尺度性,且表现出纵向分层特性和多属性场耦合的相关性。基于InterfaceGrid设计理论模型和数据结构,构建统一描述三维地质结构的InterfaceGrid数据模型,实现复杂地质结构的矢量栅格一体化表达。基于GeoAtom地原子理论构建InterfaceGrid的形式化表达框架,提出InterfaceGrid的构建流程,以及基于InterfaceGrid模型设计数据更新和空间查询算法。采用GPU光线投射和自适应采样策略实现地质网格的体可视化和在线浏览,结果显示,与SBRT模型相比较,InterfaceGrid能够更加真实地描述地质边界情况,提高了三维地质结构模型的准确性。在全球岩石圈三维网格化构建应用InterfaceGrid,验证了InterfaceGrid在大规模地质数据组织管理中的适用性,相比较PillarGrid数据量减少约1/3,更加适合基于大数据量的网络地学应用服务。

关键词

三维地质结构模型 / InterfaceGrid / 可视化 / 断层建模 / 空间查询

Key words

3D structural geological model / InterfaceGrid / visualization / fault modeling / spatial query

0 引言

前人对三维地质模型开展过许多研究。它是地质空间内地质构造和地质体对象边界的数字化表征模型,具有确定性、可视性和可修改性等特点^[1⇓⇓-4]。目前主要有结构模型和属性模型两种方法实现地质模型的构建和表达^[5-6],其中结构模型侧重于地质对象几何信息构建,主要内容是利用表达地质体的相关原始数据,在计算机中建立描述地质构造形态和地质体对象之间关系的模型^[7-8]。

三维地质结构模型目前主要用两种方法进行表达:一种是使用三角网(TIN)表达边界来划分地下空间;另一种是网格化方法,强调将地质空间离散为网格化单元,通常为常规六面体元素。其中三角网最为流行,主要原因是TIN对复杂表面的刻画能力以及在保持对象间边界拓扑一致性方面具有很好的优势,但是三角网数据规模大,影响可视化的分析的效率^[8⇓-10]。因此,复杂地质模型的网格化表达在科学研究和工程应用中发挥着越来越重要的作用,广泛用于自然资源的评估和预测、岩石物理性质的估计和地质统计学的研究^[7⇓⇓-10]。

地质模型的网格化表示方法使用不同的网格填充单元来划分空间。规则网格通过标准长方体或者平行六面体来划分空间,结构简单,在大规模数据需求下内存过大;PillarGrid^[11]的结构基于一组从地质模型的顶部到底部延伸的垂直或倾斜的柱子。柱子中的单元是由8个相邻顶点定义的六面体,顶点之间相互独立,单元可以发生形变,能够表达常见的断层,但是需要记录8个顶点,结构复杂;四面体网格^[12]自适应表达能力强,但需要记录复杂的网格信息,包括每个网格的顶点信息和连接方式,内存占用量比较大且空间查询困难;Stack-Based Representation of Terrains(SBRT)^[13]水平方向上仍是规则网格,在垂直方向上记录岩层边界位置和对应的岩性;Vertically Projected Triangulated Network(VPTN)^[9]在水平方向上做三角网剖分,垂直方向上根据地层界面切割出三棱柱,记录切割点高程和地层岩性;SBRT和VPTN的内存更为紧凑,结构简单,但是无法模拟复杂的地形现象,如Y型断层和分支断层。

随着多源异构地质数据矢量栅格一体化集成和地质数据空间查询分析需求增多,基于不规则三角网的方法不能满足要求,地质模型的网格化表达方法成为当前的研究热点^[13⇓-15]。数据采集技术的进步使地球数字化变得越来越准确,由此产生的大量数据显著影响计算资源和性能,如存储和可视化所需的内存大小,以及传输和处理。这极大地影响了整体模型的表达。同时,由于地质信息具有极强的空间非均匀性、不连续性、多尺度和纵向分层的多属性场耦合等特性,目前的网格化方法,如规则网格、PillarGrid和SBRT模型等,仍然无法实现复杂地质结构的精细表达^{[14⇓⇓⇓⇓⇓⇓-21]}。构建统一表达地质结构的一体化空间数据模型是三维地质信息科学领域需要解决的基础性问题之一。

本研究结合承担项目平台开发的空间分析和可视化的需要,提出InterfaceGrid,一种内存紧凑且利于空间分析的矢量栅格一体化三维地质结构模型表达方法。InterfaceGrid与目前提出的结构模型不仅具有同等的能力,并且从理论上能够对整个复杂地质模型的表达起到有力的补充作用。InterfaceGrid能够减少模型占用的内存空间和减少计算量,能够更准确地刻画三维复杂地质,并且易于组织管理大规模数据,整体上能够为大规模三维复杂地质数据的空间分析和可视化提供支撑。本研究提供3个应用场景,以验证InterfaceGrid模型的有效性。

1 InterfaceGrid数据模型

1.1 形式化理论表达

GeoAtom是一个元组<p,z,z(p)>,它描述一个特定位置(p)的属性值信息。一个属性和可能取值的几何被称为域和范围^[22]。p的维度不受限制,可以用于三维空间或者四维时空。GeoAtom可以组合成地理物体或者地理单元。地理数据包括点线面体,被称作地理对象,地理原子在空间离散域上聚合为地理场。

地理场的离散方式有两类:基于瓦片的地理场(TBGF)和基于采样与插值的地理场(S&IBGF)^[23-24]。在TBGF中,域中任意空间位置获得的采样值与该位置所在瓦片相关,一个瓦片中包含的所有GeoAtom都具有同样的属性值。S&IBGF使用特定的插值方法,计算场中特定点的值,例如为了获取遥感影像中点的数据,使用邻近的像素点数据双线性插值计算数据结果。

由于地理场是由一个不确定的GeoAtom集合组成的,因此有必要制定一个离散化方案,以简单表示地理场数据且易于数据管理。在地形和地质建模中,可以在地理场中采样的属性之一是特定点的地理属性,此属性可以通过标量TBGF进行映射,TBGF创建了地理原子集合,可以描述为^[22]

<p,M,f_m(p)>

式中:M={m₀,m₁,m₂,...,m_n,u}是一个简化的属性集合;u表示地理属性;p是任意空间位置;f_m是一个映射函数,定义为

f_m:R₃→M

在地形和地质建模中,存在许多离散化地理场的方法^[9⇓⇓-12]。三维规则网格是经典的地理场空间离散方案,对于每个输入点p可以明确映射到一个体素,体素的属性值即为p点的属性值;规则网格是地理场离散化的常规方法,结构简单,广泛应用于隐式地质建模中,但是可能会带来较高的内存占用。

作为一种地理场离散化方法,InterfaceGrid模型的形式化定义为InterfaceGridModel={Ω,I,R},其中I={I_i}, i=1,...,n表示模型中地质界面的数量,R={r_ij}, r_ij=<I_i,I_j>表示两个不同地质界面I_i、I_j的关系, z_i(x,y)表示地质界面I_i在x、y位置处的z值。

r_ij= \left\{\begin{array}{l}{I}_i\mathrm{和}{I}_j\mathrm{相}\mathrm{交}\iff I_i\bigcap I_j\ne \mathrm{\varnothing }\\ I_i\mathrm{位}\mathrm{于}{I}_j\mathrm{上}\mathrm{方}\iff z_i(x,y)\ge z_j(x,y)\mathrm{或}\\ I_i\mathrm{位}\mathrm{于}{I}_j\mathrm{下}\mathrm{方}\iff z_i(x,y)\le z_j(x,y)\end{array}\right.

假设I_i是单值曲面,I_i的分段函数表达为I_i=z(x,y)= \left\{\begin{array}{l}{z}_1\left(u\right),u\in {\Omega }_1\\ z_2\left(u\right),u\in {\mathrm{\Omega }}_2\\ \dots \\ z_n\left(u\right),u\in {\mathrm{\Omega }}_n\end{array}\right.,u=(x,y),定义域Ω₁∪Ω₂∪…∪Ω_n⊆Ω,其中Ω_j可以是三角形、四边形等几何图形。

如果Ω_j是规则网格,定义域Ω:m×n∈R²,空间有(m+1)×(n+1)个点,参数的节点矢量节点坐标X={x₀≤x₁≤…≤x_m},Y={y₀≤y₁≤…≤y_n},分段函数的子区间u_ij=(x_i,y_j),x∈X,y∈Y,Ω_ij=(u_i,j:u_i+1,j+1)=(x_i:x_i+1,y_i:y_i+1)。

参数曲面: \mathrm{z} ( \mathrm{x} , \mathrm{y} ) = \stackrel{m}{\sum _{i=0}} \stackrel{n}{\sum _{j=0}} z_{ij} N_i \left( x \right) N_j \left( y \right)称为一次B样条曲面。式中:z_ij称为曲面的控制点;逐次用线段连接点列z_ij中相邻两点组成的空间网格称为曲面的控制网格;N_i(x)N_j(y)是由节点矢量X,Y定义的规范B样条基函数。

本研究采用的基函数是线性插值函数。一次B样条曲面是由m×n片一次参数多项式组合而成的多项式样条曲面,曲面的控制网格S_ij=z_ijN_i(x)N_j(y)表示一个Interface,且I_i={ I n t e r f a c e_j^i}=<z_i,j,z_i+1,j,z_i+1,j+1,z_i,j+1>

推论:

(1)I_i位于I_j的上方,且I_j位于I_k的上方,则I_i位于I_k的上方等同于

z_i(x,y)≥z_j(x,y)且z_j(x,y)≥z_k(x,y),则z_i(x,y)≥z_k(x,y)

(2)I_i位于I_j的下方,且I_j位于I_k的下方,则I_i位于I_k的下方等同于

z_i(x,y)≤z_j(x,y)且z_j(x,y)≤z_k(x,y),则z_i(x,y)≤z_k(x,y)

1.2 数据结构

InterfaceGrid模型数据结构设计的基本思路是根据地质数据确定定义域,在定义域内构建全局模型,模型单元表达为单元柱(Pillar)。水平方向将网格抽象为XYGrid数组,Pillar则表达为地质界面(Interface)和空间关系(二叉树)的组合,如图1^[25]所示。

InterfaceGrid模型需要表达两部分信息:第一部分是高程和属性数据^[25-26],可以由静态数组表达;第二部分是地质界面的关系信息,直观上表现为链式二叉树。由于链式存储的内存开销较大,本研究使用静态数组表达链式二叉树,同时设计静态数组与链式二叉树之间相互转换的算法以支持数据更新。

XYGrid数组中存储网格单元存储数据位置和地质界面片段的数量。假设网格规模为M_x×M_y,XYGrid是一个M_x×M_y×2的三维数组,以行优先顺序线性化,在索引(x,y)位置,三维数组与一维数组的对应关系如下:

XYGrid[x][y][0]=XYGrid[(y×M_x+x)×2]

XYGrid[x][y][1]=XYGrid[(y×M_x+x)×2+1]

Interfaces数组存储网格单元中地质界面片段的高程数据,将浮点型高程值在整数范围内离散化,使每个高程值节省4字节的内存消耗。对于高程范围在顶和底的网格数据模型,64位操作系统中实际存储的z值计算如下:

z=(z-z_Bottom)/(z_Top-z_Bottom)×2³²-2³¹

数组Interfaces的内存大小由所有网格单元中Interface的数量决定, Interface的顺序对应于链式二叉树的先序序列。数组InterfaceProp存储Interface属性编码及右子树信息。网格索引部分的数据量取决于分辨率,地质空间中存在N个Interface时,数据部分占用N×48个字节,InterfaceGrid数据模型的空间复杂度与Interface数量呈线性关系。

2 构建流程

在InterfaceGrid数据模型中,利用基于曲面的边界表示方法,在特定空间关系约束下使用统一网格化流程和组织方式来构造、存储和管理地下曲面对象的拓扑结构和几何结构。三维地质结构模型的构建涉及支持的数据模板定义和数据预处理、空间关系抽象化和空间划分,以及构造规则制定等方面^[21]。基于InterfaceGrid理论表达将模型的构建划分为以下具体步骤,构建流程图如图2所示。

(1)制定建模参数。计算层面控制点数据的AABB包围框作为定义域Ω,根据区域地质背景和建模源数据分析结果设置水平方向分辨率R_xy,三维地质属性场离散化为网格单元Grid,称为基网格,如图3所示。

(2)空间关系定义。使用已有的钻孔、地调、物探资料和经验推断出地层层序关系和断层的接触关系,将层面对象抽象为节点并以二叉树森林定义空间关系。

(3)地质界面网格化。在地质勘查数据中获得地下地质构造的采样点信息,可以利用采样点信息了解地质构造和地质对象的几何形态,但是地质构造中其他点的信息是未知的^[22]。因此,需要利用地质勘查获取的采样点信息拟合网格化曲面来预测未知点信息。本研究中采用薄板样条插值算法利用控制点数据拟合地质界面,已知某个地质界面控制点集合为{(x₀,y₀,z₀),(x₁,y₁,z₁),…,(x_n,y_n,z_n)},则薄板样条插值函数为^[27-28]:

\mathrm{\phi } ( \mathrm{x} , \mathrm{y} ) = \stackrel{n}{\sum _{j=1}} A_j {d_j}^2 \mathrm{l} \mathrm{o} \mathrm{g} d_j + a + b x + c y = z

式中:a,b,c为系数;A∈R^1*n;d_j为待插值点(x,y)与地质界面控制点中第j个点的距离;插值点的采样密度与InterfaceGrid模型的水平方向分辨率一致。

(4)一体化建模。步骤(3)中生成的网格化地质界面投影到地质网格中形成Interface,空间关系包括地层之间、断层之间、地层与断层之间的关系,由(2)定义的空间关系森林逐个网格点嵌入基网格生成Interface,建立整体模型的过程。整体模型在水平方向上为均匀网格,纵向上记录建模工区内的地质界面,地层之间的层序关系、尖灭关系,地层与断层以及断层之间的接触关系由(2)中定义的二叉树表示。

3 关键算法

3.1 地质界面数据更新

考虑界面之间复杂的交切关系,在添加和删除地质界面数据时,需要由静态数组构建二叉树,执行常规二叉树数据更新算法实现。本研究设计基于InterfaceGrid数据结构的二叉树创建算法(BuildBiTree)。

BuildBiTree函数的功能是使用给定的先序遍历序列和各先序结点对应的右子树序列构造链式二叉树。其中,先序遍历序列按照根节点、左子树和右子树的顺序排列,右子树序列存储节点右子树在先序遍历序列中的位置。函数伪代码如下:

BuildBiTree函数的复杂度可以表达为O(n),其中n为树中结点的数量。该函数通过递归的方式实现了对树的构造。在每一次递归中,函数都会创建一个新的结点,并记录节点属性。接着,函数会根据rightTreeOrder数组的信息递归地构造结点的右子树,再递归地构造结点的左子树。由于在每一次递归中都会对Flag数组进行标记,保证了每个结点只会被构造一次。因此,对于n个结点的树,该函数的时间复杂度为O(n)。因此,该函数的时间复杂度与二叉树的节点个数呈线性关系。表1所示为本研究二叉树创建方法与先序和中序序列构建二叉树方法的复杂度对比分析。

3.2 空间查询算法

在可视化过程中,需要快速有效地获取到模型空间任意位置或者任意区域的属性信息。根据InterfaceGrid数据结构特点,可以设计一种稳健高效的空间数据查询方法QueryPoint。

QueryPoint总体上是一个递归函数,首先确定待查询点point所在单元柱pillar。根据索引数据获得地质界面数据,即Interfaces和InterfaceProp。假设待搜索单元柱的地质界面(Interface)数量为N, index表示正在处理的interface索引值,初始值为0,最大值为N-1;invalidIndex记录无效查询区域,即在先序序列中,invalidIndex索引位置以外的Interface数据属于已被排除的范围;top和bottom分别记录待查询点的地质界面索引位置。

对于正在处理的地质界面(curInterface)数据,使用双线性插值计算查询点位置在Interface上对应的高度值(z)。比较待查询点与Interface的空间位置关系,若待查询点在较上面的位置,更新bottom为curInterface,index向右滑动即数值加1。同时由于当前处理的Interface已经作为底部,该Interface空间位置以下的所有Interface将无须判断,因此更新invalidIndex为位于curInterface下面的地质界面索引(curInterface的右子树位置);否则更新top为curInterface,index更新为位于curInterface的右子树索引。

递归结束条件有两个,满足其一则程序结束:(1)index与invalidIndex指示位置重合;(2)index为无效值(空间位置超出模型数据底部边界)。最终,待查询点位于top和bottom对应的两个地质界面组成的地层中,该地层的属性信息即查询点的属性信息。算法伪代码如下:

单元柱空间可以看作是二叉树空间分区,二叉树空间分区是一种常见的空间划分方法,它递归地将n维空间划分为关于划分超平面的子空间^[14]。Interface是单元柱中用于划分空间的平面。通常,在二叉树空间分区中需要建立链式二叉树才能进行查询,从数据构建一棵平衡二叉树的时间复杂度通常是O(nlogn)。本研究提出的查询算法的显著优势在于无须构建二叉树,利用两个序列数组进行递归查询。QueryPoint函数的空间复杂度主要是由递归调用所引起的系统栈空间消耗。每次递归调用都会将当前的函数栈帧压入系统栈中,并在递归结束后将该函数栈帧弹出,因此系统栈的最大深度取决于平衡二叉树的高度h,时间复杂度为O(log(n)),其中n是平衡二叉树中的节点数。表2所示为本研究查询方法与常规构建二叉树查询方法的复杂度对比分析。

4 建模实例及效果分析

4.1 模型对比测试

三维地质构造的表达能力是数据模型可视化效果的重要评价标准,SBRT模型由于其数据量和空间查询的优势被广泛应用于在线可视化领域^[13⇓-15]。本研究应用测试数据三维网格模型对比验证SBRT模型和InterfaceGrid模型的体绘制效果。它共有6个岩层,岩性属性在三维网格内呈现层状分布,存在岩层尖灭现象。水平方向上划分网格规模为64×64,网格间距为10 m×10 m。测试环境为:CPU为Intel Core i7, 显卡为Nvidia Geofore MX150。采用谷歌浏览器,体绘制效果如图4所示。

由图4可见,SBRT模型作为规则网格的纵向压缩,继承了规则网格的特点,体绘制时在岩层边界处会产生明显的锯齿或台阶问题,影响可视化效果,不能准确描述地质情况;InterfaceGrid模型直接记录地质界面在当前网格点位置的深度值且单元内使用矢量表达,模型岩性边界无锯齿问题,放缩时模型形态一致,在一定分辨率上准确还原地质真实情况,更符合真实数据效果,提高了三维地质结构模型表达的准确性。在绘制效率方面,SBRT模型的绘制效率是40帧率,InterfaceGrid数据模型的绘制效率是25帧率。由于InterfaceGrid模型的数据量比SBRT模型大,InterfaceGrid模型体绘制的显示帧率比SBRT模型低。随着GPU硬件的升级,计算机显存提升明显,因此内存导致的一定的帧率降低是可接受的^[11⇓-13]。

4.2 全球球面网格剖分测试

测试数据LITHO1.0模型的球面是一个正二十面体球面三角网,在每个三角形单元的径向上描述了岩石圈的9个亚层(水、冰、上层沉积物、中层沉积物、下层沉积物、上地壳、中部地壳、下地壳和岩石圈地幔) ^[29]。针对全球岩石圈模型多源异构地质空间信息整合和地质大数据一体化组织需求,本研究应用来自LITHO1.0模型的9个亚层数据,利用相同分辨率的Citcoms网格^[30-31]建立了PillarGrid模型和InterfaceGrid模型,验证InterfaceGrid模型。PillarGrid和InterfaceGrid均为柱状网格,基于一组从地质模型的顶部到底部延伸的垂直或倾斜的柱子。PillarGrid的一个单元由其8个相邻顶点定义,与InterfaceGrid相比,它在纵向上重复记录了相邻单元的共用顶点,两种不同结构的单元对比如图5所示。图6为LITHO1.0模型的全景岩石圈层显示效果,两种模型外观无明显差异。两种模型的单元数量和内存空间占用量如表3所示。

InterfaceGrid模型中的单元数量比PillarGrid模型的单元数量增加了约8%,这是由于InterfaceGrid的单元为四边形面,PillarGrid模型的单元为六面体,同一位置纵向上面单元与体单元表达范围存在差异。InterfaceGrid模型中的内存占用量比PillarGrid模型减少了约36.8%,对于不同区域的模型,这一比例可能会有所不同,与建模区域的地层条件有关。一般来说,在相同的运行环境中,平均层数越大,InterfaceGrid模型与PillarGrid模型的比值越小,总体上InterfaceGrid模型中的内存占用量约是PillarGrid模型的2/3。

5 结论

根据地质空间数据矢量栅格一体化集成以及对数据查询、分析和计算的需求,采用网格化的空间数据表达方法已成为地学数据集成的主要研究趋势。针对规则网格、角点网格和广义三棱柱等无法实现复杂地质结构的表达问题,本研究考虑地质结构和属性的地下分布具有极强的非均匀性、不连续性和空间多尺度性,且表现出纵向分层特性和多属性场耦合的相关性,提出了InterfaceGrid数据模型,应用它统一描述地质结构和多属性场的地下分布,以实现地质模型矢量栅格的一体化表达,从而实现地质空间数据查询、分析、可视化及应用。本研究测试验证表明,与SBRT模型相比,InterfaceGrid模型的可视化结果边界无锯齿问题,能够更加真实地描述地质边界情况,提高了三维地质结构模型的准确性和可靠性。

本研究使用Citcoms全球球面网格剖分方法构建地球表面到莫霍面9个岩石圈分层的InterfaceGrid模型数据并实现可视化。实验显示,InterfaceGrid模型可以实现全球化圈层结构的网格化描述,同时InterfaceGrid模型相比较PillarGrid数据量更小,更加适合基于网络的大数据地学应用服务。

综上,作为矢量栅格一体化模型,InterfaceGrid兼具栅格模型和矢量模型的特点和优势: 一方面,InterfaceGrid结构简单,空间水平坐标隐含在网格索引中,纵向上以二叉树序列构建,无额外内存产生且查询效率高,利于进行空间分析计算和在线可视化;另一方面,InterfaceGrid以空间四边形表达地质界面,可视化时无锯齿问题,一定程度上减少了网格尺寸导致的模型与真实数据之间的空间误差,从而提高了模型的准确性和可靠性,有利于大尺度精细建模。

后续研究将采用自适应多尺度网格和并行处理提升InterfaceGrid的存储管理和计算的性能,也可以将基网格由规则网格推广到三角形网格,进一步丰富InterfaceGrid的理论框架。

参考文献

[1] 毛先成, 张显洋, 刘占坤, 胶西北金矿集区地球化学多重分形模型与异常识别[J]. 物探化探计算技术, 2022, 44(2): 213-223.

[2] 袁峰, 张明明, 李晓晖, 成矿预测: 从二维到三维[J]. 岩石学报, 2019, 35(12): 3863-3874.

[3] 周永章, 左仁广, 刘刚, 数学地球科学跨越发展的十年: 大数据、 人工智能算法正在改变地质学[J]. 矿物岩石地球化学通报, 2021, 40(3): 556-573, 777.

[4] 潘懋, 方裕, 屈红刚. 三维地质建模若干基本问题探讨[J]. 地理与地理信息科学, 2007, 23(3): 1-5.

[5] MALLET J L. Geomodelling[M]. New York: Oxford University Press, 2002.

[6] WANG Z G, QU H G, WU Z X, et al. Formal representation of 3D structural geological models[J]. Computers and Geosciences, 2016, 90(PA): 10-23.

[7] 曾佐勋, 樊光明. 构造地质学[M]. 3版. 武汉: 中国地质大学出版社, 2008.

[8] 唐骥, 蒋潇, 姜雪莲, 矿体三维可视化建模技术在成矿模式分析中的应用[J]. 地质科技通报, 2023, 42(5): 273-284.

[9] ZHOU C Y, DU Z C, OUYANG J W, et al. A 3D geological model and cutting algorithm based on a vertically projected triangulated network[J]. Computers and Geosciences, 2020, 143: 104562.

[10] CAUMON G, COLLON-DROUAILLET P, LE CARLIER DE VESLUD C, et al. Surface-based 3D modeling of geological structures[J]. Mathematical Geosciences, 2009, 41(8): 927-945.

[11] PEYROT J L, DUVAL L, PAYAN F, et al. HexaShrink, an exact scalable framework for hexahedral meshes with attributes and discontinuities: multiresolution rendering and storage of geoscience models[J]. Computational Geosciences, 2019, 23(4): 723-743.

[12] JACKSON M D, HAMPSON G J, SAUNDERS J H, et al. Surface-based reservoir modelling for flow simulation[J]. Geological Society, London, Special Publications, 2014, 387(1): 271-292.

[13] BENES B, FORSBACH R. Layered data representation for visual simulation of terrain erosion[C]// Proceedings spring conference on computer graphics. Budmerice, Slovakia: IEEE, 2001: 80-86.

[14] HOULDING S W. 3D geoscience modeling: computer techniques for geological characterization[M]. Berlin, Heidelberg: Springer Science and Business Media, 2012.

[15] 李青元, 张丽云, 魏占营, 三维地质建模软件发展现状及问题探讨[J]. 地质学刊, 2013, 37(4): 554-561.

[16] LI F Y, GAO C L, LIU Y Q, et al. Integrated multi-scale reservoir data representation and indexing for reservoir data management and characterization[J]. Computers and Geosciences, 2020, 138: 104433.

[17] HILL K M, SWIFT J N, HOWARD C, et al. 3D visualization of subsurface objects from La Brea tar pits, Los Angeles, CA[J]. Digital Applications in Archaeology and Cultural Heritage, 2021, 20: e00167.

[18] 蔡思敏, 任伟中, 冯亮, 基于GTP-TEN的复杂地质体三维混合建模[J]. 岩石力学与工程学报, 2023, 42(2): 441-449.

[19] SUN J X, LENZ D, YU H F, et al. MFA-DVR: direct volume rendering of mfa models[J]. Journal of Visualization, 2024, 27(1): 109-126.

[20] QIN R F, FENG B, XU Z N, et al. Web-based 3D visualization framework for time-varying and large-volume oceanic forecasting data using open-source technologies[J]. Environmental Modelling and Software, 2021, 135: 104908.

[21] YOUNG G, KRISHNAMURTHY A. GPU-accelerated generation and rendering of multi-level voxel representations of solid models[J]. Computers and Graphics, 2018, 75: 11-24.

[22] 熊祖强, 贺怀建, 夏艳华. 工程地质三维建模及分析系统设计研究[J]. 岩石力学与工程学报, 2007, 26(增刊2): 4176-4182.

[23] 陈良, 高成敏. 快速离散化双线性插值算法[J]. 计算机工程与设计, 2007, 28(15): 3787-3790.

[24] GRACIANO A, RUEDA A J, FEITO F R. A formal framework for the representation of stack-based terrains[J]. International Journal of Geographical Information Science, 2018, 32(10): 1999-2022.

[25] PASYANOS M E, MASTERS T G, LASKE G, et al. LITHO_1.0: an updated crust and lithospheric model of the Earth[J]. Journal of Geophysical Research: Solid Earth, 2014, 119(3): 2153-2173.

[26] GRACIANO A, RUEDA A J, FEITO F R. Real-time visualization of 3D terrains and subsurface geological structures[J]. Advances in Engineering Software, 2018, 115: 314-326.

[27] 孙黎明, 刘禹杉, 张睿卓, 基于几何细分的三维地质模型自适应精细化构建方法[J]. 岩土工程学报, 2023, 45(增刊1): 244-248.

[28] 牛露佳, 王双威, 曾义文, 基于规则网格的复杂断层网络处理与建模[J]. 地质论评, 2023, 69(5): 1980-1990.

[29] YU J Q, WU L X, ZI G J, et al. SDOG-based multi-scale 3D modeling and visualization on global lithosphere[J]. Science China: Earth Sciences, 2012, 55(6): 1012-1020.

[30] ZHONG S J, MCNAMARA A, TAN E, et al. A benchmark study on mantle convection in a 3-D spherical shell using CitcomS[J]. Geochemistry, Geophysics, Geosystems, 2008, 9(10): 1-32.

[31] 陈建平, 李靖, 谢帅, 中国地质大数据研究现状[J]. 地质学刊, 2017, 41(3): 353-366.

基金资助

广东省重点研发计划项目(2020B1111370001)

国家重点研发计划项目(2022YFF0801201)

国家自然科学基金重点项目(U1911202)